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CHAPTER  1 


INTRODUCTION 


1.1  The  Generalized  Minimum  Cost  Network  Flow  Problem 


In  recent  years  there  has  been  an  increasing  awareness  that 
network  flow  programming  is  a  powerful  modeling  and  problem  solving 
technique  that  can  be  applied  to  a  wide  range  of  physical  and 
conceptual  situations.  The  network  flow  programming  model 
encompasses  a  variety  of  problem  classes,  including  the  shortest 
path  problem,  the  maximum  flow  problem,  the  pure  minimal  cost  flow 
problem,  the  convex  network  flow  problem,  and  the  generalized 
minimum  cost  flow  problem.  Among  these  classes,  the  generalized 
network  is  the  problem  that  is  considered  in  this  research.  The 
generalized  transportation  problem  and  generalized  assignment 
problem  are  special  cases  of  a  generalized  minimum  cost  flow 
problem. 

The  network  of  a  generalized  minimum  cost  flow  problem  is  a 
directed  graph  defined  by  a  set  of  arcs,  M,  with  ordered  pairs  of 
nodes  (i,j)  as  elements  indexed  by  k.  For  each  arc  in  the 
generalized  minimum  cost  flow  problem  there  is:  a  cost  h^  for  each 
unit  of  flow  which  passes  through  the  arc;  a  lower  bound  c^,  which 
is  the  minimal  allowable  flow  that  an  arc  can  carry  ;  a  capacity  cl 


which  is  the  maximum  amount  of  flow  that  an  arc  can  carry  ;  and  a 
gain  factor  a^  which  determines  if  flow  is  to  be  lost,  generated,  or 
conserved  on  the  arc.  Each  node  in  the  minimum  cost  flow  problem  is 
either  a  supply  node  where  flow  enters  the  network,  a  demand  node 
where  flow  leaves  the  network,  or  a  transshipment  node  where  no  flow 
enters  or  leaves  the  network.  The  required  quantities  of  flow 
entering  or  leaving  the  network  at  the  nodes  are  called  external 
flows.  A  positive  external  flow  enters  the  network  at  a  node,  and 
a  negative  external  flow  leaves  the  network  at  that  node. 

The  goal  of  a  generalized  minimum  cost  flow  problem  is  to 
determine  how  a  commodity  should  be  delivered  through  the  arcs  of  a 
network,  that  is  to  choose  the  arc  flows  ffc,  so  that  shipment  cost 
will  be  minimized  and  two  kinds  of  constraints  are  satisfied.  The 
two  constraints  to  be  satisfied  are  :  1)  that  flow  be  conserved  at 
each  node,  which  is  to  say  that  flow  leaving  a  node  from  the  arcs  of 
a  network  minus  flow  entering  the  node  on  the  arcs  of  the  network 
must  equal  the  external  flow  at  the  node;  and  2)  that  the  flow  on 
each  arc  be  between  its  lower  bound  and  capacity.  The  algebraic 
representation  of  the  model  can  be  expressed  as  a  bounded  linear 
programming  problem  with  a  constraint  for  each  node  and  a  variable 
for  each  arc. 

The  algebraic  model  is  stated  as  follows: 


S.T 


i*  1  )  •  •  •  jtl“l 


i  fk  i  Ck 

or  in  matrix  form 


k  £  M 


Min  HF 
S.T.  AF  -  b 

C  <  F  <  C 


where 

f.  is  the  amount  of  flow  on  arc  k 
k 

hfc  is  the  shipping  cost  for  each  unit  of  flow  on 
arc  k 

c^  is  the  minimum  amount  of  flow  on  arc  k 

c.  is  the  maximum  amount  of  flow  on  arc  k 
k 

a^  is  the  gain  parameter  on  arc  k 

is  the  external  flow  at  node  i  (positive  for 
incoming  flows,  negative  for  outgoing 
flows) 


MQi  is  the  set -of  arcs  which  originate  at  node  i 


is  the  set  of  arcs  which  terminate  at  node  i 
F,H,C,C,  and  b  are  respectively  the  collections 
of  f^,  h^,  c^,  c^,  and  b^.  A  is  the  node 
arc  incidence  matrix  of  the  generalized 
minimum  cost  flow  model.  An  example  of  the 
A  matrix  will  be  shown  in  Chapter  4. 


1.2  Parametric  Programming 

An  Important  part  in  the  analysis  of  any  large  mathematical 
programming  problem  is  the  study  of  how  the  solution  changes  with 
variations  of  the  problem  data.  There  may  be  questions  concerning 
the  accuracy  of  the  data  or  a  simple  desire  to  see  how  the  solution 
changes  as  elements  of  the  solution  change.  This  problem  gives  rise 
to  the  area  of  sensitivity  analysis.  When  performing  sensitivity 
analysis,  one  determines  the  region  over  which  the  restrictions  may 
by  changed  so  as  to  maintain  an  optimal  solution.  This  analysis  may 
be  extended  to  parametric  programming,  where  a  series  of  optimal 
solutions  are  generated  as  the  conditions  are  continuously  changed. 
This  allows  the  analyst  to  be  able  to  answer  a  series  of  "what  if" 
questions  about  the  problem.  The  flexibility  created  by  allowing 
changes  in  any  part  of  the  problem  data  leads  to  a  much  more 


powerful  analytic  technique. 

In  particular  the  analyst  may  wish  to  vary  some  of  the 
right-hand  side  restrictions.  These  might  Include  variations  of  the 
external  flows  (b)  or  the  capacity  restrictions  (C).  We  can  state 
the  parametric  programming  problem  for  varing  the  external  flows  as 
follows: 


Min  HF 

S.T.  AF  -  b  +  0  b 
C  <  F  <  C 
0  <  Q  <  9  max 

where  0  is  the  change  to  be  accomplished  in  the  external  flow  vector 
b.  Similarly  we  can  state  the  parametric  programming  problem  for 
the  variation  of  the  capacity  vector  as 


Min  HF 
S.T.  AF  •  b 

A 

C  <  F  <  C  +  X  C 


0  <  X  <  X 


max 


where  X  is  the  change  to  be  accomplished  in  the  capacity  vector  C 
We  can  also  give  a  concise  verbal  description  of  the  parametric 
programming  problem  as  follows: 


6 


Given:  An  opcimal  solution  to  a  generalized 

minimum  cost  network  flow  problem  and  a  set 
of  parameters  of  the  problem  which  are  to 
be  varied. 

Problem:  Determine  a  series  of  optimal 

solutions  to  the  original  problem  as  the 
parameters  are  continuously  varied. 

By  Choice  of:  How  much  flow  to  send  along  each 
arc  of  the  network. 

Subject  to:  External  flow  restrictions  must  be 
met;  conservation  of  flow  at  each  node  must 
be  maintained;  lower  and  upper  bounds  on 
each  arc  must  not  be  violated. 


The  significance  of  this  solution  technique  is  that  it  will 
allow  the  analyst  to  examine  a  range  of  optimal  solutions  for  the 
problem  being  modeled.  Secondly  it  opens  the  avenues  to  the 
solution  of  problems  which  have  Imbedded  networks  as  subproblems  in 
which  either  capacities  or  external  flows  are  being  varied.  The 
solution  algorithm  could  take  the  previous  iterations'  optimal 
solution  and  iterate  by  use  of  the  parametric  sensitivity  algorithm 
which  will  be  presented. 


1.3  Study  Outline 

In  this  research,  an  attempt  is  made  to  develop  a  parametric 
programming  algorithm  which  applies  to  the  generalized  network  flow 
problem.  A  detailed  review  of  past  literature  relative  to  the 
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history  of  linear  programming,  generalized  network  programming,  and 
the  solution  technique  is  presented  in  Chapter  2.  In  Chapter  3  a 
summary  of  linear  programming  theory  relevant  to  this  problem  is 
discussed.  Chapter  4  is  a  review  of  pertinent  network  flow 
terminology  and  theory.  Chapter  5  presents  a  detailed  description 
and  explanation  of  the  solution  technique  as  it  is  applied  to  the 
generalized  network  flow  problem.  Finally,  in  Chapter  6,  two  new 
dual-incremental  flow  algorithms  are  developed.  Computational 
results  comparing  these  codes  to  two  primal  codes  and  a  published 
dual-incremental  code  are  presented. 


CHAPTER  2 


LITERATURE  REVIEW 


2.1  Introduce ion 


The  generalized  minimum  cost  network  flow  problem  is  an 
important  class  of  network  flow  problems  in  which  efficient  solution 
algorithms  exist.  The  generalized  network  flow  problem  is  a  special 
case  of  a  bounded  variable  linear  program.  Since  linear  programming 
has  played  such  an  Important  role  in  the  development  of  network  flow 
programming,  we  first  review  the  history  of  linear  programming 
before  further  discussion  of  generalized  network  flow  literature. 

2,2  Linear  Programming  Background 

Some  of  the  original  works  of  linear  programming  theory  can 
be  traced  back  to  1823  when  Fourier  experimented  with  homogeneous 
linear  systems.  These  methods  were  similar  to  the  later  discovered 
simplex  algorithm  (28).  During  this  same  period  the  works  of  Gauss 
and  Jordan  also  established  some  basic  techniques  for  solving  linear 
systems  of  equations.  P.  Gordan  appears  to  be  the  first  to  recognize 
that  a  linear  system  of  homogeneous  equations  has  an  associated  dual 
homogeneous  problem  and  that  non-negative  variables  in  these  systems 
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of  equations  possesses  a  solution  with  at  least  one  variable 
positive  if  the  dual  possesses  no  solution  with  strict  Inequalities. 
Steimke  expanded  this  result  to  include  all  positive  variables  (21). 
Applying  Fourier's  principles,  Minkowski  concluded  that  in  a 
homogeneous  system,  the  general  system  can  be  formed  as  a 
non-negative  linear  combination  of  finite  extreme  point  solutions. 
Vllle,  'Von  Neuman,  and  Morgensteln  presented  relationships  between 
primal  and  dual  solutions  in  the  late  1930's  and  early  1940's,  but 
it  was  Kantorovich,  in  1939,  who  came  extremely  close  to  developing 
the  first  linear  programming  algorithm.  He  outlined  solution  methods 
to  a  class  of  problems  consisting  of  simultaneous  outputs  (49).  His 
initial  approach  was  based  on  an  initial  feasible  dual  solution. 

The  most  influential  work  appeared  to  be  that  of  Neyman  and 
Pearson  (63),  in  1936,  whose  development  of  some  duality  and  linear 
inequality  theories  led  to  the  significant  discovery  of  Dantzig 
(21).  The  first  primal  algorithmic  procedure  for  linear  programming 
problems  was  the  simplex  method  developed  by  Dantzig  in  1947. 
Shortly,  thereafter.  Von  Neuman  expressed  the  additional  importance 
of  duality  in  solving  linear  systems  (78).  In  1951  Tucker,  Kuhn, 
and  Gale  developed  a  duality  theorem.  Also  in  1951,  Dantzig  and 
Orden  proved  that  the  optimal  solution  for  the  primal  provided  a 
corresponding  optimal  for  the  dual.  Lemke,  1954,  discovered  the  dual 
simplex  algorithm.  He  also  added  an  important  note  on  the 


complementary  slackness  conditions  between  the  primal  and  dual 
solutions  (56).  Charnes,  et  al.  (18,19)  pioneered  the  application 
of  linear  programming  techniques  to  the  industrial  world,  besides 
publishing  numerous  articles  which  contributed  to  linear  programming 
theory.  Continued  search  for  new  algorithmic  approaches  resulted  in 
Orchard-Hay's  composite  simplex  which  is  used  in  those  cases  where  a 
solution  was  found  to  be  neither  primal  or  dual  feasible  (64). 

These  results  formed  the  basis  of  linear  programming. 

During  the  1950-1970's  with  the  rise  of  improved  computer 
technology,  much  of  the  research  done  was  to  improve  the 
computational  efficiency  of  the  simplex  on  large  scale  programming 
applications.  The  introduction  of  the  revised-simplex  technique,  of 
Dantzlg,  Orchard-Hays  and  others  led  to  more  efficient  numerical 
methods.  Dantzig  and  Van  Slyke  (22)  improved  the  computational 
efficiency  of  a  special  class  of  linear  programs  with  generalized 
upper  bounds  (GUB)  constraints.  Hellerman  and  Rarlck  (41)  and 
Forrest  and  Tomlin  (27)  worked  on  numerical  stabilty  for  reinverting 
the  basis.  Schrage  (68,69)  and  Todd  (75)  developed  special 
algorithms  for  applications  of  the  simplex  method  to  linear 
programs  with  variable  upper  bounds  (VUB). 

The  refinement  of  data  storage  and  retrieval,  as  well  as 
more  efficient  computational  methods,  have  been  a  springboard  in 
recent  years  for  the  methods  mentioned  above  to  be  applied  to  a  wide 


field  of  complex  large  scale  linear  programs. 

2.3  Generalized  Networks 

The  first  computational  procedures  for  the  generalized 
network  problem  related  to  the  generalized  transportation  problem. 
This  application  has  a  bipartite  structure  and  uncapacitated  arcs. 
Dantzlg  (21)  and  Charnes  and  Cooper  (17)  described  early  attempts  at 
adapting  the  simplex  method  to  this  special  problem.  Eisemann  (24), 
Lourle  (37),  and  Balas  and  Ivanescu  (11)  lmplefnented  the  use  of 
primal  procedures  in  relation  to  the  stepping  stone  method  of  the 
pure  transportation  problem.  Eisemann  (24)  additionally  discussed 
the  use  of  a  dual  solution  approach  to  the  problem.  All  these 
methods  use  the  matrix  representation  of  the  transportation  problem 
and  are  based  on  the  bipartite  nature  of  the  graph. 

Jewell  (47)  Introduced  a  procedure  similar  in  concept  to  the 
out  of  kilter  algorithm  (26)  to  solve  the  generalized  minimum  cost 
flow  problem  on  a  general  network  with  capacities.  Minleka  (61) 
modified  Jewell's  algorithm  to  guarantee  finite  termination. 

Johnson  (48)  suggested  the  use  of  the  three  pointer 
representation  of  the  basis  in  the  generalized  network  problem. 
Maurras  (59)  and  Glover,  Kllngman,  and  Stutz  (37)  utilized  a  three 
pointer  representation  to  Implement  primal  algorithms  for  the 


generalized  minimum  cost  flow  problem.  Additional  computational 
simplifications  were  made  to  the  generalized  transportation  problem 
by  Glover  and  Klingman  (36).  Jensen  and  Bhaumik  (44)  described  a 
dual  incremental  approach. 

Hultz  and  Klingman  (43)  provided  an  improved  dual  feasible 
start  and  made  computational  comparisons  on  pure  network  test 
problems.  Glover,  et  al.,  (38)  presented  an  application  paper  on  the 
use  of  generalized  networks  and  computational  testing  with  their 
generalized  network  code  NETG.  Elam,  Glover,  and  Klingman  (29) 
presented  a  primal  simplex  algorithm  which  is  specially  designed  to 
reduce  the  number  of  degenerate  pivots.  Jensen  and  Barnes  (46) 
developed  a  generalized  primal  code,  called  PGAINS,  and  a  dual 
incremental  code  call  INCREMG.  Adolphson  (2)  extented  the  preorder 
thread  to  generalized  models  in  his  computer  code  AGENNET.  Glbby, 
et  al.,  (32)  explored  different  pivoting  strategies  to  improve  the 
efficiency  of  primal  codes. 

The  interest  in  the  area  of  applications  of  generalized 
networks  continues  to  expand,  as  does  the  size  of  the  problems  which 
are  being  solved  by  these  solution  methods.  We  have  merely  surveyed 
the  major  accomplishments  in  the  area  of  generalized  networks. 


2.4  Sensitivity  ifoalysis  and  Parametric  Programming 

During  the  history  of  linear  programming  many  economists  and 
specialists  have  met  with  failure  when  they  have  attempted  to 
introduce  linear  programming  into  their  operations.  This  has  been 
based  frequently  on  one  of  the  following  factors: 

1.  The  uncertainty  and  inaccuracy  of  the  initial  data. 

2.  The  problems  of  evaluation  and  interpretation  as  well  as 
application  and  exploitation  of  the  results  in  practice. 

The  problems  mentioned  in  1.  gave  rise  to  the  area  of  sensitivity 
analysis  in  which  one  may  test  in  what  regions  the  right-hand  side 
restrictions  may  be  changed  to  maintain  the  optimality  of  the 
solution.  The  earliest  use  of  parametric  programming  as  applied  to 
linear  programs  dates  back  to  Gass  and  Saaty  (30,31).  Parametric 
programming  enables  us  to  compute  all  existing  optimal  basic 
solutions  in  relation  to  their  dependence  on  the  values  of  the 
components  of  the  right-hand  side.  Gal  (29)  lists  over  500 
references  in  his  bibliography  in  relation  to  either  parametric  or 
sensitivity  analysis  of  linear  programs. 

The  same  problem  of  uncertainty  and  inaccuracy  of  the 


initial  data  should  lead  to  similar  problems  in  network  solution 
analysis  by  today's  managers.  However  these  problems  have  not  led 
to  a  great  deal  of  research  interest  in  this  area.  Thus,  while  there 


exist  a  large  number  of  references  in  the  linear  programming  area, 
there  are  relatively  few  references  in  the  area  of  parametric 
analysis  for  network  flow  problems.  It  appears  that  the  initial 
research  in  a  related  are  was  done  by  Srinivasan  and  Thompson 
(72,73)  when  they  applied  operator  theory  to  parametric  programming 
for  the  transportation  problem.  In  a  series  of  articles  Balachandran 
and  Thompson  (5, 6, 7, 8)  apply  operator  theory  to  the  generalized 
transportation  problem.  Srinivasan  and  Thompson  (74)  develop  a 
theory  of  cost  operator  algorithms  for  transportion  problems.  In 
these  methods  parametric  programming  is  use  to  solve  the 
transportation  problem  to  optimality.  Balachandran  (10)  developes 
an  operator  theoretic  approach  for  a  generalized  transportation 
problem  with  stochastic  demands.  Balachandran,  Srinivasan  and 
Thompson  (9)  demonstrate  the  application  of  iheir  operator  theory  of 
parametric  programming  for  both  transportation  and  generalized 
transportation  problems.  Ahuja  (3)  developes  a  ranging  analysis  for 
pure  capacitated  transshipment  problems.  Adolphson  (1),  in  a 
private  communication  to  the  author,  describes  a  procedure  for 
ranging  analysis  of  the  right-hand  side,  cost,  and  gains  factors  for 
a  generalized  network  flow  problem  which  has  a  computation  time  of 
order  m  (where  m  is  the  number  of  arcs). 


CHAPTER  3 


D 

- 

I 


BACKGROUND  IN  LINEAR  PROGRAMMING 

3.1  Introduction 

Since  a  capacitated  generalized  network  is  a  special  case  of 
a  linear  program,  it  seems  appropriate  to  discuss  the  fundamental 
ideas  behind  the  linear  programming  simplex  method.  In  this  chapter 
we  will  first  review  the  bounded  variable  simplex  method,  give  a 
brief  discussion  of  the  dual-simplex  method,  and  lastly  develop  the 
background  for  perturbation  of  the  right-hand  side  of  the 
constraints. 

3.2  Statement  of  the  Problem 


Using  matrix  notation,  the  linear  programming  problem  with 
simple  upper  bounds  may  be  stated  in  the  form: 

Min  Z  -  HTF 
S.T.  AF  -  b 

0  <  F  <  C. 


(3.1) 

(3.2) 

(3.3) 
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While  the  above  notation  is  not  standard  for  linear  programming,  it 
follows  the  network  notation  used  by  Ford  and  Fulkerson  (26).  F  is 
the  vector  of  decision  variables.  H  is  the  vector  of  variable 
costs.  The  right-hand  side  values  are  represented  by  the  vector  b, 
and  C  is  the  vector  of  capacities  or  upper  bounds.  There  are  n 
variables  and  m  equations . 

Equation  (3.1)  is  called  the  objective  function  and 
equations  (3.2)  and  (3.3)  are  called  constraints.  A  is  an  m  by  n 
matrix  with  n  >  m.  We  will  assume  the  equations  are  linearly 
Independent,  that  is,  we  do  not  have  any  equation  which  is  a 
multiple  of  one  or  more  of  the  other  equations.  Inequality 
constraints  may  be  made  into  equality  constaints  as  in  equation 
(3.2)  by  use  of  "slack  variables"  which  are  added  to  the  left-hand 
side  and  whose  values  represent  the  difference  between  the  original 
left-hand  side  and  the  right-hand  side.  They  are  called  "logical 
variables";  the  original  variables  are  called  "structural  variables" 
as  they  have  more  physical  significance. 

The  main  concept  behind  upper  bounding  is  that  any  nonbasic 

variable  fR  may  either  be  at  zero  of  its  upper  bound  c^.  This 
eliminates  the  need  to  represent  the  upper  bound  on  a  variable  as  a 
separate  constraint,  thus  saving  a  row  in  the  constraint  matrix  for 
each  upper  bound.  This  reduces  m,  the  size  of  the  basis  matrix, 
with  the  resultant  savings  in  both  running  time  and  storage 


requirements.  The  price  paid  for  such  savings  is  that  the  steps  of 
the  simplex  algorithm  become  slightly  more  complicated. 

We  can  use  equation  (3.2)  to  solve  for  m  variables  called 
the  basic  variables  in  terms  of  the  other  (n-m)  variables  called  the 
nonbasic  variables.  We  partition  A  into  the  submatrices  B,  , 


*  ’  1  8  I  *l»  I  *ob  1  <3-‘> 

where  B  is  an  m  by  m  nonsingular  submatrix  formed  from  the  columns 
of  the  basic  variables  called  the  basis  matrix;  R^g  and  R,,g 
constitute  the  columns  of  the  nonbasic  variables  at  zero  and  the 
columns  of  the  nonbasic  variables  at  their  upper  bounds, 
respectively. 

We  can  express  equation  (3.2)  in  the  form 

BFB  +  RLBFLB  +  RUBFUB  *  b  (3.5) 


where  F  has  been  partioned  into 
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F  - 


(3.6) 


corresponding  Co  Che  parcicioning  of  A. 

EquaCion  (3.5)  can  be  rearranged  Co  Che  form 


"b  ‘  1  •  Vu  '  Wn  (3-7> 

and  since  B  is  nonsingular 

Fb  .  »-‘b  -  B-\A,  -  (3.8) 

We  have  chus  solved  for  m  of  Che  variables  (F  )  in  Cerms  of  Che 

O 

ocher  (n-m)  variables  (F^,  Fyg),  where 

Fg  is  called  Che  vecCor  of  basic  variables 

Ff„  is  called  Che  vecCor  of  nonbasic  variables 
uo 

ac  zero 

F  is  called  Che  vecCor  of  nonbasic  variables 
aC  Che  upper  bound. 

A  basic  soluCion  is  defined  as  one  in  which  Che  nonbasic  variables 
are  sec  aC  eicher  cheir  upper  bound,  C^,  or  zero.  A  basic  feasible 
soluCion  is  one  in  which  all  Che  cerms  of  che  veccor  Ffi  are 


nonnegacive 


3.3  The  Transformed  Problem 


As  we  partitioned  F  into  (Fg,  FLfi,  Fyg) »  we  can  also 
partition  H  into  (Hg,  H^g,  Hyg).  We  can  now  rewrite  equation  (3.1) 
as 


Min  Z 


HB  FB  +  HLB  flb  +  HUB  fub* 


(3.9) 


From  equation  (3.8)  we  have  an  expression  for  F„  in  terms  of  FT _  am 

o  LB 

Fyg  and  can  write  the  objective  function  as 


Z 


+  H 


<»'lb  - 

LB  FLB  +  HUB  FUB* 


-  b‘1Wbb> 


(3.10) 


Rearranging  terms 

Z  •  »BTB_1b  -  (H„TB-\b  -  HlbXb 

-  <hbTb'\b  *  "obT)fob-  <3-"> 

Thus  we  can  express  our  original  problem,  equations  (3.1  -  3.3),  in 
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z  -  “b1''1"  -  (hbTb'\b  -  "lbXb 

(Hg  B  Ryg  -  Hyg  )Fyg  (3.12) 

S.T.  Fb  -  B”ib  -  B’1^  -  B^RygFyg  (3.13) 

FB  >  °- 

This  Is  called  Che  transformed  problem,  or  sometimes  referred  to  as 
the  current  tableau.  If  we  let 

6  -  B-1b 

and  ‘ 

T  -1 

ir  -  HLB 

a 

the  transformed  problem  can  be  expressed  as 

Min  2  .  KbTS- 

^ 71  rub  *  kub  ^fub 

S.T.  F,  -  B  -  B-\bFu  -  B-\,Fub 

F,>  o- 

Since  F^g  ■  0,  Fg  takes  on  the  numerical  value  of  6  -  B  an<* 

the  objective  function  takes  on  the  value  of  HgS  -  (ff  R^g  -  »UB)FUB- 
We  leave  the  right-hand  terms  of  the  equations  in  the  transformed 
problem  even  though  F^g  ■  0,  since  they  indicate  what  happens  to  Z 


(3.14) 

(3.15) 


and  F_  as  one  of  the  elements  of  F  Is  Increased  from  zero  or  one 
a  LB 

of  the  elements  of  F^  Is  decreased  from  its  upper  bound. 


3.4  Conditions  for  Optimality 


The  row  vectors 

(3.16) 

(3.17) 


(  *T*u  -  "ibT> 
(  *\b  -  hubT) 


which  premultiply  FLB  and  F^g  In  equation  (3.14),  the  objective  row 
of  the  transformed  problem,  are  called  the  vectors  of  reduced  costs, 
since  they  indicate  how  much  the  objective  Z  decreases  as  F 

LB 

increases  from  zero  or  the  amount  Z  increases  as  Fyfi  decreases  from 
CUB* 

Let  us  denote  the  jth  element  of  the  reduced  cost  vectors, 
equations  (3.16-3.17)  by 


djLB  "  (zjLB  ‘  hjLB^ 
djUB  “  (zjUB  “  hjUB^ 


for  nonbasic  variables  at  zero 
for  nonbasic  variables  at  upper 


bounds 
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where  h  is  the  jth  element  of  H  and  h  is  the  jth  element  of 

J  LB  *  LJ3  j  UO 

hub*  and 


2jLB-  <*\b>)  ■  ”TajLB  (3'l8> 

ZjOB  "  ^UB^j  "  2  ajOB 


where  a.TO  and  a,fTll  are  the  jth  columns  of  R.  „  and  R„_,  the  nonbaslc 

j  LB  jUd  LB  Uo 

portion  of  A. 

We  see  that  if  for  f^LB  the  corresponding  reduced  cost  is 
positive,  Z  will  decrease  upon  increasing  f...  above  zero,  while  if 

J  LB 

for  the  corresponding  reduced  cost  is  negative,  Z  will  decrease 

on  decreasing  f  below  c. . 

jun  j 

Thus  if  our  transformed  problem  represents  a  basic  feasible 
solution  with  Fb  -  6  -B~ \gFyg  >  0,  Fyj  -  0,  and  FyB  -  CUB,  then 
the  further  condition  for  optimality  is  that  all  elements  of  the 
reduced  cost  vectors  be  negative  (or  zero)  for  F  and  positive  (or 
zero)  for  F  . 

Summarizing  these  conditions: 

Feasibility  S i  “  b"1riuBF1UB  -  0  i-l,....« 

dj  ^  ° 
dj  >° 


Optimality 


jt  LB 
je  UB 
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where 

LB  -  {je  R  j  fj  is  at  0} 

UB  -  {je  R  |  fj  Is  at  Cj}  . 


3.5  Pivoting 


We  consider  the  possibility  of  increasing  only  one  nonbasic 
variable  from  zero  or  decreasing  from  its  upper  bound,  while  all 
others  stay  at  zero  or  their  upper  bounds.  If  more  than  one  reduced 
cost  is  positive  for  a  nonbasic  variable  at  zero  or  negative  for  a 
nonbasic  variable  at  its  upper  bound,  we  have  an  arbitrary  choice. 

It  is  customary  in  most  introductory  texts  to  consider  the  nonbasic 
variable  with  the  "largest"  reduced  cost  since  this  will  decrease 
the  value  of  the  objective  function  the  most  per  unit  change.  There 
are  more  elaborate  ways  to  pick  the  entering  variable,  but  we  will 
not  discuss  them  here  (see  reference  62). 

Writing  out  the  transformed  problem  column  by  column,  it 
takes  the  form 


Mlnz 

^n-m^FR^n-m 

S.T.  Fg  -S  “  ••• 


...  -  d  (F_) 
q  R  q 


-<*  (F„)  '  •« 

q  R  q 


(3.20) 
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where 


an-m^FR^n-m 


a. 


B_1a. 


j  -  1, 


tn-m  • 


(3.21) 


(3.22) 


Here  we  have  chosen  a  nonbaslc  variable  (F_)  to  Increase  from  zero 

R  q 

or  to  decrease  from  Its  upper  bound.  Since  some  of  the  nonbaslc 
variables  may  be  at  their  upper  bound,  S  does  not  represent  the 
current  value  of  F  .  Rather,  we  have 

D 


F  »  B  -  I 
B  j  eUB 


aj(Vj 


(3.23) 


where  UB  is  the  set  of  nonbaslc  variables  at  their  upper  bounds. 
To  dlstinquish  between  current  and  new  values,  let  (F  ) 

O  1 

denote  the  current  value  of  the  ith  basic  variable,  and  F^  denote 

its  new  value  as  (F  )  changes. 

R  q 

As  (fR)q  Increases  (decreases),  each  element  In  F^  will 

change : 

1.  If  a.  Is  positive,  F.  will  decrease 

lq  1 

(increase) 

2.  If  Is  negative,  F^  will  increase 
(decrease) . 
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We  need  Co  consider  three  possible  events: 


a)  A  basic  variable  may  reach  its  lower  bound 
(zero)  first 

b)  A  basic  variable  may  reach  its  upper  bound 
first 

c)  A  nonbasic  variable  may  reach  the  other  end 
of  its  range  before  (a)  or  (b)  occur. 


Suppose  (F  )  is  increased  from  zero. 
R  q 

For  case  (a)  we  have 


or 


<Vi  -%<FR>q  >  ° 


(F  )  <  min  {(F  )  /  a  )  . 

Rq“*i|cL>0  Bi  lq 


For  case  (b)  we  have 


or 


(FVi  ~  “iq^R^q  -  Ci 


<Vq<  “1"  «Ci  -  <Vi>/  "V 

i|  aiq<° 


For  case  (c)  we  have 


(3.24) 

(3.25) 


(3.26) 

(3.27) 
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<FR>q  <  Cq.  (3.28) 

Thus  (F  )  must  take  on  the  smallest  value  of  (3.25) ,(3.27) ,  or 
R  q 

(3.28).  Equivalently,  we  can  write 


(F  )  -min  {  min  {(FR). 

Rq  i|a.  >0  Bi  iq 

iq 

min  {(C  -  (F  )  )/  -a.a},  C  } . 

i|a.  <0  1  81  lq  q 

iq 


(3.29) 


We  oust  also  consider 
The  current  values  of 


the  case  where  (F  )  is  decreased  from  C  . 

R  q  q 

the  basis  variables  are  given  by 


B  i 


j£UB, 


(3.30) 


Thus  for  any  new  value  of  (F  )  ,  we  have 

r  q 


<Vi  -  <'»>,>• 


(3.31) 


If  we  define  the  decrease  in  the  nonbasic  variable  by  <5  »  - 

(F  )  ,  then  for  case  (a)  we  have 

R  q 


^B^  +aiq(5)  >  0 


(3.32) 
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6  <  min  { ( F_) .  /-a  }  .  (3.33) 

81  ^ 


For  case  (b)  we  have 


(3.34) 

(3.35) 

For  case  (c)  we  have 


or 


6  <  min  {(C.  -  (F  )  )  /  ot  }  . 
-1|%>0  1  81  “> 


6  <  C  . 

—  q 


(3.36) 


Combining  the  results  of  equations  (3.33),  (3.35),  and  (3.36)  we 
have 


6  «  min  {  rain  {(F_).  / 

*1  V° 

min  {(C  -  (F  )  )/a  } ,  C  }  (3.37) 

i  I  a  >0  q  q 

iq 


and  (F„)  ■  C  -6  .  We  should  note  if  (c)  occurs  no  change  of  the 

R  q  q 

basis  is  required.  However,  we  need  to  record  the  fact  that  (F  ) 

K  q 


A 
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has  gone  to  tie  other  end  of  its  range.  The  new  basic  variables  F 

B 

written  in  terms  of  the  old  basic  variables  F  are 

B 


✓ 


an(+C  ) 

q  —  q 


and  the  appropriate  sign  is: 


+  if  (F  )  increased  from  zero 
k  q 

-  if  (F_)  descreased  from  C  . 
R  q  q 


If  (a)  or  (b)  occurs,  for  example  if  i  ■  p,  F  is  adjusted  according 

a 

to 


/ 


<fb>i  'W, 

(FR^q 


i*p;  (F_)  increased  from  zero 
k  q 

ifp;  ( F_)  decreased  from  C 
k  q  q 

i-p. 


At  this  point  a  non-basic  variable  has  a  positive  value  and  a  basic 

variable  has  a  value  of  zero  or  is  at  its  upper  bound.  Recall  that 

in  the  partitioned  equation  (3.7),  (F  )  corresponds  to  the  pth 

R  p 

column  of  B.  The  pivot  operation  is  completed  when  the  tableau  is 

rearranged  so  (F„)  is  basic  and  (F„)  is  non-basic.  This  means 
R  q  Bp 

redefining  the  partition  of  A  so  a^  replaces  the  pth  column  of  B 
corresponding  to  entering  the  basis  and  (F^)p  leaving  the 

basis.  Note  the  fact  that  a  record  must  be  kept  of  which  nonbasic 
variables  are  at  their  upper  bounds. 
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3.6  Steps  of  the  Simplex  Method 

The  procedures  outlined  in  sections  3.2  through  3.5  form  the 
building  blocks  of  the  bounded  variable  simplex  method.  While  there 
exist  numerous  variants  of  this  method,  most  computer  codes  use  the 
revised  simplex  method.  Thus  in  all  large  scale  applications  all 
the  elements  of  -B  *R  of  the  transformed  problem  are  never  stored  or 
computed.  We  only  need  two  columns  of  the  transformed  problem  (  8  ■ 

B  *  b  and  a  «  B  *a  ).  We  assume  we  have  a  basic  feasible  solution 
with  some  representation  of  B  1  and  the  current  right-hand  side  6  * 
B  1b.  The  steps  of  the  simplex  are  as  follows: 

STEP  1:  Produce  a  pricing  vector 

ttT-HbTB-1.  (3.38a) 

STEP  2:  Price  out  the  nonbasic  columns 

T 

a)  Evaluate  z,  ■  it  a,  (3.38b) 

b)  Evaluate  d^  ■  1  h^.  (3.38c) 

STEP  3:  Select  the  entering  nonbasic  variable  by  choosing 


a 
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a)  The  most  positive  d^  among  the  nonbasic 
variables  at  zero 

b)  The  most  negative  among  the  nonbasic 
variables  at  upper  bounds 

c)  If  none  exists,  stop;  otherwise  choose  q 
corresponding  to  the  reduced  cost  of  the 
largest  magnitude  of  a)  or  b) . 

STEP  4:  Update  the  entering  column  by  evaluating 

a  -  B_la  .  (3.38d) 

q  q 

STEP  3:  Find  the  leaving  basic  variable  by  selecting 

(^r)_  -  min  {  min  {(Fg),  ^  a  ia*  • 

"  *1  ct  lq>0  ‘  1 

,wLn  <(Ci  -  (FB>i>/  -aio>.  ca>  ( 3 . 38e ) 

or 

6-  min  {  min  {(F^  /  -at  } , 

il  aiq<0  q 

min  {{C±  -  (FgM/a  CQ}.  (3.38f) 

^V0 
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.  a)  If  none  exists,  stop  with  an  unbounded 
solution 

b)  Otherwise  let  the  minimum  correspond  to  i  * 
P* 


STEP  6:  Pivot  by  adjusting  so 

(F„)  increased 
x  q 

from  zero 

(F0)  decreased 
k  q 

from  Cq  (3.38g) 

of  B. 

Return  to  step  1. 


<Vi  -  aiq<FR>q  i+p; 


(FB)i'aiq(6) 


<Vq 


i“P 


and  aq  replaces  the  pth  column 


3.7  The  Dual-Simplex  Method 


This  procedure  applies  the  simplex  method  to  the  dual 
problem  while  working  with  the  primal  tableau.  Suppose  we  have  the 
standard  primal  tableau  expressed  in  the  form 


Min  Z  -  H„T  B-  dl(FR),  -  ...  -  dq(Fa)q  -  ... 
^n-m^FR^n-m 
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S.T. 


?B  “  ®  “al^*Vl  ~  **•  _Qlq^FR^q  ”  *** 


“  “n-m^R^n-m* 


(3.39) 


For  ease  of  exposition  this  explanation  describes  the 
dual-simplex  method  for  a  linear  program  without  bounds.  We  assume 
the  tableau  is  dual  feasible  (optimal)  but  not  necessarily  primal 
feasible.  That  is,  all  the  reduced  costs  are  less  than  or  equal  to 
zero  (dj  <_  0) ,  but  one  or  more  of  the  elements  of  the  vector  b  is 
negative.  Then  the  steps  of  the  dual-simplex  method  are: 


STEP  1:  Choose  the  pivot  row  p  by  selecting 


6  *  min  {  6 . }  .  (3.40) 

P  i  |  ^<0  1 


a)  If  none  exists,  stop  with  optimal  and 
feasible  solution 

b)  Otherwise  let  the  minimum  correspond  to 
l*p,  say  (F_)  is  the  leaving  basic 
variable.  p 


STEP  2:  Choose  the  entering  nonbasic  variable. 


We  must  note  that  the  reduced  cost  of  (F_)  in  the  new 

B  q 

tableau  if  (FR)^  replaces  (FR)^  in  the  basis  is  given  by 
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dj  -  -dj  /  (aj)p  .  (3.41) 

The  ratio  test  now  becomes 

0  *  min  {-d.  /  (a.)  }  (3.42) 

jl  ^ j)p<0  J  3 

for  j  -  q.  Noting  the  fact  that  if  any  higher  ratio  were 

used,  say  for  j»r  and  (F_)  entered  the  basis,  and  if  (-d  / 

k.  r  r 

<Vp>  >  (“dq  '  (Vp>’  thSn 

(Vp  1  dr  <  0  (3*34) 

and  optimality  would  be  lost.  Thus  when  (F„)  enters  the 

k  q 

basis 

dj  "  dj  +  ®  foj)p  0.44) 

STEP  3:  Check  for  infeasibility 

If  there  is  no  element  (<*j)p  <  0  for  any  of  the 

nonbasic  variables  (FR)^,  then  no  nonnegative 

value  of  (F„).  can  cause  (F_)  to  be  feasible; 

K  J  ®  P 


Therefore,  stop.  The  problem  Is  infeasible. 
Otherwise  go  to  step  4. 


STEP  4:  Perform  the  pivot  operation 


As  in  the  primal  method,  exchange  (F_)  with 

k  q 

(F_)  in  the  basis  and  return  to  step  1. 

B  q 


There  are  two  major  motivations  for  using  the  dual-simplex 


method. 


a)  Adding  extra  rows  to  the  primal  problem. 
Using  the  dual-simplex  method  this  is 
achieved  with  essentially  the  same  ease  as 
adding  an  extra  column  to  the  primal  and 
pricing  out  the  new  column  with  the  current 
basis. 

b)  It  is  also  used  in  parametric  programming 
of  the  right-hand  side  vector. 


3.8  Perturbation  of  the  Right-Hand  Side  Vector  C 


During  the  development  of  parametric  programming,  most  of 


the  research  effort  has  been  spent  on  analyzing  the  problem 
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Min  Z  - 

htf 

/\ 

(3.45) 

S.T. 

AF  -  b  +0  b 

(3.46) 

in  the  range 

0  <  0  <  0 

—  max 

(3.47) 

F  >  0. 

(3.48) 

However  little  analysis  has  been  spent 

on  the  development  of  the 

parametric  analysis 

for  the  problem 

- 

Min  Z  -  HTF 

(3.49) 

S.T.  AF  - 

b 

A 

(3.50) 

C  < 

F  <  C  +  X  C 

(3.51) 

°  <  x  <  xma„ 

max 

where  we  parametrically  vary  the  upper  bounds  of  the  variable  F, 
with  x  being  a  scalar. 

A  capacitated  generalized  network  is  a  special  case  of 

problem  (3.49).  Therefore,  in  order  to  do  parametric  analysis  on 

the  capacitated  generalized  network,  we  must  first  analyze 

parametric  changes  in  the  implicit  right-hand  side  constraints. 

A 

Thus  suppose  the  vector  C  is  replaced  byC+^C,  ^^0. 

/\ 

This  means  the  capacity  constraints  are  changed  along  the  vector  C. 
Since  the  right-hand  side  of  the  primal  problem  is  the  objective  of 
the  dual  problem,  changing  the  implicit  capacity  constraints  could 


be  analyzed  as  changing  the  objective  function  of  the  dual  problem. 
However,  we  shall  now  discuss  the  effects  of  changing  the  capacity 
constraints  directly  in  the  primal  problem. 

Recalling  from  section  3.2,  we  can  partition  matrix  A  into 
(B,R).  We  denote  f  corresponding  to  the  ith  column  of  B  by 
and  let  r' be  the  set  of  indices  of  the  columns  of  R.  Let 

Qt  j  ■  B_1aj  jeR' 

and 

6  -  B-1b  . 

We  can  write  equation  (3.50)  in  solved  form  as 


. . 


(3.52) 


If  we  define 

LB  -  {j  e  r'I fj  is  at  Cj} 

UB  ■  { j  £  RX I  f  J  is  at  Cj  +  x  c j  } 

Z  -  {j  e  r"|  f  is  at  0} 


then  we  can  write  (3.52)  as 
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Bi 


B'1  !  b  +  Z  <*(-<•)  +  Z  a  (_x£ ‘ 
j  e  LB  J  J  jeUBJ  J 


-  c, 


If  we  let 


-I  , 

P  -  B  [  -  Z  a  .c. 
°  jeUB  3  3 


and  lec 


--1 


a  -  b  [  b  -  Z  a  c.  -I  a  c.  ] 
°  j  e  LB  3—3  je  UB  3  3 


then  we  can  rewrite  (3.53)  as 


PBi-0lo*Cl,X*  J2Vj  ‘-1 . . 


Now  we  consider  varying  X  to  satisfy 


-l-°io+PloX  i  cl  +  X  ci 


where  LB^,  UB^  are  the  bounds  on  the  basic  variable  Fg^. 
From  (3.54)  we  require  F  to  satisfy 


°  <.  (  aiQ  -  ct)  +  (  p^)  X  1-1,..., m 


) 

(3.53) 


(3.54) 


(3.55) 


r  * 

(  alo  -  cA)  +  (plo  “  ci)X£0  i»l . . 


(3.56) 


In  general  we  consider  increasing  X  from  zero  to  some  upper  limit 
✓\ 

X.  We  define 


Q,  -  {It  (aio  -  £l)  >0,  Plo  <  0! 

Q2  *  Ul  (  °io  '  cl>  i  °-  <  plo  "  V  >  °>  • 


(3.57) 


As  X-is  increased,  (3.55)  may  be  violated  if  and  only  if  ieQj  and 
(3.56)  may  be  violated  if  and  only  if  i£  Qj.  Thus  let 


X1  "  ?enQl{(aio  "~i)  '  ~P*°} 

A 

X2  -  min  {-(a  lo  -  ct)  /  (  p lo  -  ct)} 


X  ■  min  (  X  j,X  p. 


(3.58) 


If  the  minimum  is  attained  in  (3.58)  forX^,  let  p  ■  p^.  Thus  Fgp 
will  leave  the  basis  at  its  lower  bound  (k»l)  or  its  upper  bound 
(k»2).  This  leads  to  two  cases. 

CASE  1:  X  »Xj.  In  this  case  Fgp  leaves  the  basis  by  going  to  its 
lower  bound.  If  we  were  to  increase  X  further,  it  must  be  true  from 


L 


(3.52),  we  must  either 


(a)  increase  a  nonbasic  variable  which  has 
j  <  0  and  j  e  LB  \j  Z 

or 

(b)  decrease  a  nonbasic  variable  which  has 
(Xjj  >  0  and  j  e  UB  U  Z. 

If  such  a  variable  can  be  found,  say  F_  ,  perform  a  simplex  pivot 

Kq 

and  make  FR^  basic  while  Fg  enters  the  set  of  nonbasic  variables  at 
their  lower  bound. 

CASE  2:  X  ■  In  this  case  Fg  leaves  the  basis  by  going  to  its 

upper  bound.  If  we  increase  X  further,  it  also  must  be  true  from 
(3.52)  that  we  must  either 

(a)  Increase  a  nonbasic  variable  which  has 
a  p j  >  0  and  j  t  LB  (J  Z 

or 

(b)  decrease  a  nonbasic  variable  which  has 
otpj  <  0  and  j  e  UB  U  Z. 

Again  if  such  a  variable  can  be  found,  say  F_  ,  make  the  appropriate 

Kq 

simplex  pivot  by  letting  F„  become  basic  while  F_  enters  the  set 

Kq  op 

of  nonbasic  variables  at  their  upper  bound. 
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If  no  element  B  r  of  the  appropriate  sign  can  be  found  in 

pq 

the  "blocking  row"  p  at  any  step,  the  algorithm  terminates. 

An  algorithm  to  accomplish  the  above  procedure  would  have 
the  following  general  steps: 


STEP  0:  Solve  the  original  problem  to  optimality. 
Set  i  ■  I,X  *  0,  and  X^  ■  0. 


STEP  1:  Determine  the  element  B  r  of  the 

pq 

appropriate  sign.  If  none  can  be  found,  stop; 
the  current  basis  is  optimal  for  all  values  of  X 
>  Xj^  or  any  further  iterations  will  cause  the 
problem  to  become  infeasible.  Otherwise 
calculate  X  . 


STEP  2:  Let  Xj_^«-  Xj^  and  Xj_«-  X.  For  Xe 

[  X^j.  X  1  the  current  basis  will  remain 

optimal.  Remove  F  from  the  basis  and  pivot 
ap 

F  into  the  basis.  Update  the  elements  of  the 
Kq 

simplex  tableau  and  the  objective  function. 
Return  to  step  1. 


1 
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3.9  Perturbation  of  the  Right-Hand  Side  Vector  b 

If  in  equation  (3.50)  we  replace  b  by  b  +X  b,  we  note  the 

term  H^b  -  HR  will  not  be  affected,  that  is  dual  feasibility  will 

not  be  affected.  The  only  change  is  that  B~*b  will  be  replaced  by 
-1  ^ 

B  (b  +  X  b)  and  accordingly  the  objective  function  becomes 

-1  ^ 

HfiB  (b  +X  b). 

As  long  as  HgB^Cb  +  X  b)  remains  nonnegative,  the  current  basis  will 
remain  optimal.  We  can,  therefore,  determine  the  value  of  X  at 
which  another  basis  becomes  optimal.  Let  S  ■  {il  B”*b^  <  0}.  If  S 

■  <!>,  then  the  current  basis  is  optimal  for  all  values  X  ^  0. 

Otherwise,  as  before  let 

X  -  min  {B  V  /  -B  b  )  . 
i  G  S  1 

A 

Let  X  j  »  X  .  ForXc[  0,  X^  ],  the  current  basis  is  optimal,  where 

-1  ~ 

FB  -  B  *(b  +X  b). 

-1  ^ 

At  X  j  the  right-hand  side  is  replaced  by  B  (b  +  X  b) ,  F^  is 
removed  from  the  basis,  and  an  appropriate  variable,  according  to 


the  dual-simplex  criteria,  enters  the  basis.  The  process  is 
repeated  to  find  the  new  range  [  X  ^,X  ^  1  over  which  the  new  basis 

A 

is  optimal  t  X  ■  X^  ].  We  terminate  the  process  when  S  is  empty,  in 
which  case  the  current  basis  is  optimal  for  all  values  of  X  greater 
than  or  equal  to  the  last  value  of  X,  or  else  when  all  the  entries 
in  the  row  whose  right-hand  side  dropped  to  zero  are  nonnegative. 

In  this  case  no  feasible  solutions  exist  for  all  values  X  greater 
than  the  current  value. 


CHAPTER  4 


BACKGROUND  IN  NETWORK  FLOW  PROGRAMMING 


4.1  Introduction 


To  better  understand  the  solution  method  to  be  developed  in 
Chapter  5,  we  shall  first  discuss  some  of  the  concepts  for  the 
solution  of  a  generalized  minimum  cost  flow  problem.  A  detailed 
explanation  of  the  concepts  discussed  below  can  be  found  in  the  text 
by  Jensen  and  Barnes  (46). 

To  Illustrate  a  typical  generalized  minimum  cost  flow 
problem,  consider  the  network  model  depicted  in  Figure  4.1.  The 
algebraic  representation  of  the  model  is 

min  H  F 

S.T.  A  F  -  b 

F  <  C 

F  >  0 


where 


44 


B 


(Flow. Capacity. cost,  gain) 
[Fixed  external  flow] 


rk  with  Optimal  Solution 
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H  -  (2,20,1,12,2) 
FT  -  (3, 1,0, 1,0. 5) 


BT  -  (4, 0,0,-. 375) 

CT  -  (3, 4, 1.5, 1,1. 2). 


4.2  Network  Structure 


The  structure  of  a  network  model  Is  defined  by  nodes  and 
arcs.  A  node  1  Is  an  element  of  the  set  of  nodes,  N  »  [1,2,.  . 
.,1,.  .  . ,nj.  An  arc  may  be  defined  by  an  ordered  pair  of  nodes 
(l,j)  or  as  an  element,  say  k,  of  the  set  of  arcs  M  »  [1,2,.  .  . ,k 
.  . ,m] .  Thus  an  arc  may  be  Identified  as  an  arc  k,  arc  k(l,j),  or 
(i,j).  An  arc  k(i,j)  Is  said  to  originate  at  node  1  and  terminate 
at  node  j.  Alternatively,  1  Is  called  the  origin  node  and  j  Is 
called  the  terminal  node  of  arc  k.  We  Identify  the  origin  and 


terminal  lists 


where  and  t^  are  Che  origin  and  terminal  nodes,  respectively,  of 
arc  k.  The  collection  of  nodes  and  arcs  is  a  directed  network 
written  D  *  [N,M] ,  where  the  values  n  and  m  and  the  vectors  0  and  T 
are  sufficient  to  completely  define  the  network. 

A  variable  quantity  that  characterizes  most  of  the  problems 
to  be  considered  is  arc  flow.  Written  as  f^»  arc  flow  usually 
models  some  physical  quantity  such  as  the  flow  of  fluid,  flow  of 
people,  or  flow  of  money.  Let  f^  be  the  flow  in  arc  k  at  its  origin 
and  ffc  is  the  flow  in  the  arc  at  its  terminal  node.  In  a 
generalized  network,  flow  is  not  conserved  in  the  arc.  With  a 
nonzero  gain  parameter,  a^, 

fk  *  V  k 


If  a^  *  l  the  flow  is  conserved;  if  a^  <  1,  flow  is  lost;  and  if  a^ 
>  1  flow  is  gained  in  the  arc.  The  flow  vector  is  defined  as 
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where  Fl  indicates  the  transposed  F,  that  is  a  column  vector. 

Associated  with  each  flow  in  an  arc  is  a  cost.  The  arc  cost 
is  a  function  of  the  flow  and  is  written 


W  -  Vk 


where  h^  is  an  arc  parameter  giving  the  arc  cost  per  unit  flow.  The 
cost  vector  is  a  row  vector,  written  as 


H  *  [  hj,  hj ,...| h^  ] . 


The  arc  cost  is  a  function  of  only  the  flow  in  the  arc  and  not  a 
function  of  flow  in  the  other  arcs.  The  network  cost  is  the  sum  of 
the  arc  costs.  Thus 


HF  -  l  Vk‘ 
k-1  *  K 


The  arc  capacity  c^  is  an  upper  bound  for  the  flow  on  each 
arc.  The  arc  capacity  implies  the  constraint: 


fk  i  ck* 


The  arc  capacity  vector  is  a  column  vector,  written  as 


T 

C  I  ci*  c2» »cm  1 • 

In  many  applications  a  common  requirement  is  that  the  arc 
flow  be  at  least  as  great  as  some  lower  bound.  We  therefore  provide 
the  lower  bound  c^  as  an  input  parameter  for  each  arc  to  constrain 
the  flow  ffc  such  that 


The  lower  bound  vector  is  also  a  column  vector  and  is  written  as 

CT  -  <  £2’-“’£«  J* 

A  gain  is  the  parameter  that  models  a  linear  increase  or 
decrease  in  flow  as  it  passes  through  an  arc.  If  a^  *  1,  we  have  no 
change  in  the  flow  as  it  passes  through  the  arc.  If  •  1,  for  all 
k,  we  have  a  pure  network  problem.  When  0  <  a^  <  1,  flow  is 
decreased  as  it  passes  through  the  arc.  When  a^  >  l,  flow  is 

l 

increased.  It  is  possible  to  have  negative  gains,  but  the  physical 
interpretation  is  not  obvious.  Zero  gains  are  not  allowed. 
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External  flows  represent  connections  to  the  world  outside 
the  system  being  modeled.  There  are  two  kinds  of  external  flows  for 
a  node:  the  fixed  flow  and  the  slack  flow.  The  fixed  flow  at  node 
i,  b^,  enters  the  network  if  b^  >  0  and  leaves  the  network  if  b^  < 

0.  Any  feasible  solution  stipulates  network  flows  such  that  all 
fixed  external  flows  are  satisfied.  The  slack  external  flow  is  a 
variable  quantity  which  ranges  from  zero  to  some  upper  bound  defined 
for  the  node.  The  value  of  the  slack  external  flow  is  to  be 
determined  by  the  solution  procedure.  An  external  slack  input  is 
represented  by  an  arc  to  Che  node  from  a  specially  designated  node 
called  the  slack  node.  An  external  slack  output  is  represented  by  an 
arc  from  the  node  to  the  slack  node.  The  slack  node  is  different 
from  all  other  nodes  in  that  the  slack  node  has  the  capability  to 
absorb  or  generate  any  required  quantity  of  flow. 

A  feasible  flow  F  conserves  flow  at  all  nodes  of  the  network 
except  the  slack  node.  With  slack  arcs  defined  as  above,  all  flows 
except  fixed  external  flows  are  arc  flows.  Conservation  of  flow  for 
each  node  implies  that  the  total  arc  flow  leaving  a  node  minus  the 
total  arc  flow  entering  the  node  must  equal  the  fixed  external  flow 
at  that  node.  Algebraically,  the  conservation-of-f low  constraint 
for  a  general  node  i  is 


It  is  possible  to  express  any  network  flow  problem  with 
nonzero  arc  lower  bounds  as  an  equivalent  problem  with  all  zero 
lower  bounds.  Therefore,  we  delete  the  lower  bound  from  the 
parameter  list  by  making  the  following  transformation  for  each  arc 
k(i,j)  with  nonzero  lower  bound  c^  : 


a) 

b) 

c) 

d) 

e) 


replace  c^  by  zero 

replace  c^  by  c£  •  cfc 

replace  ffe  by  f£  -  ffc 

replace  bt  by  b^  -  b^ 

replace  b^  by  b^  *  bj 


“£k 


-  C. 


“£k 

+  W 


This  transformation  is  performed  one  arc  at  a  time  for  all  arcs  with 
nonzero  lower  bounds.  Each  new  transformation  uses  the  updated 
parameter  values  that  resulted  from  the  previous  step. 
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4.4  Expanded  Networks 


Ic  is  convenient  from  the  viewpoint  of  obtaining  optimal 
solutions  to  the  generalized  network  flow  programming  problems  to 
define  an  expanded  network  which  may  be  obtained  directly  from  the 
original  network. 

We  define  the  mirror  arc,  -k,  for  each  ke  M  such  that  arc  -k 

connects  the  same  nodes  as  arc  k,  but  has  opposite  direction. 

o(-k)  -  t(k) 

t(-k)  -  o(k). 

The  expanded  network,  Dg  -  [  N, 

[  N,M  ],  but  its  arc  set  contains  not  only  the  arcs  of  D,  but  also 

all  of  the  mirror  arcs  as  well.  That  is,  if  M  •  [1,2 . m]  then  Mg 

«(l,2,...,m,-l ,-2, . . . , -m ] . 

If  the  expanded  network  parameters  are  denoted  with 
asterisks,  their  general  definition  may  be  written  as 

* 


Mg  ] ,  has  the  same  node  set  as  D  ■ 


Forward  arcs :  h. 


*  h. 


Mirror  arcs:  h_^  *  -h^  /  a^ 

a-k  “  l/  V 

Flows  in  Che  expanded  network  imply  changes  in  Che  flows  of  Che 
original  network  by  Che  relation 

*  ,  *  ,  *  * 

A  k  “  fk  "  (f-k  )(a-k  }  “  f  ‘  f-k  '  \ 

where  is  the  change  in  flow  in  the  original  network.  The 
procedures  will  assure  that  a  forward  and  mirror  arc  will  never 
simultaneously  have  nonzero  flow. 


4.5  The  Marginal  Network 

*  * 

The  marginal  network  D  -  [  N,M  ]  has  the  same  node  set  as 
* 

D  and  D£.  The  arc  set  0  consists  of  that  subset  of 
admissible  arcs.  A  forward  arc  is  admissible  if  the  flow  is  not 
already  equal  to  the  arc  capacity,  that  is  if  the  flow  can  still  be 
Increased  on  that  arc  in  the  original  network.  A  mirror  arc  is 
admissible  if  the  flow  on  the  corresponding  forward  arc  is  not  equal 
to  zero;  that  is  if  the  flow  can  be  decreased  on  the  associated  arc 
in  the  original  network. 


Mg  that  are 


4.6  Basis  and  Basic  Solutions 


The  basis  for  a  generalized  minimum  cost  flow  problem  will 
consist  of  n-1  arcs  chosen  so  that  the  associated  n-1  columns  are 
linearly  independent.  For  a  pure  network  problem,  such  a  selection 
forms  a  directed  spanning  tree.  For  the  generalized  problem,  a 
selection  of  n-1  linearly  independent  columns  may  result  in  a  more 
general  form.  In  particular,  a  basis  for  a  problem  may  contain  one 
or  more  cycles. 

To  illustrate  these  concepts  consider  the  network -of  Figure 
4.1.  In  this  figure  no  slack  node  has  been  explicitly  defined.  For 
this  example  we  arbitrarily  choose  node  1  as  the  slack  node.  The 
constraint  matrix  formed  by  the  conservation-of-flow  equations  is 


-al 

0 

0 


ARCS 
2  3 

0  1 

-a„  -a. 


4 

1 

0 


5 

0 

1 


0  ”a4  ~a5 


2  NODES 

3 

4  . 


Note  that  the  constraint  associated  with  the  slack  node  has  been 
deleted.  This  constraint  is  redundant  and  must  be  deleted  to  obtain 
a  set  of  independent  rows.  Our  problem  is  to  choose  a  set  of  three 
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independent  columns  from  this  matrix.  The  set  will  be  independent 
if  the  square  matrix  formed  by  the  columns  has  a  nonzero 
determinant.  For  example,  if  we  choose  the  columns  whose  associated 
arcs  form  a  tree,  arcs  1,  3,  and  4,  the  basis  matrix  becomes 


The  determinant  of  the  basis  matrix  is 


|  b|  -  "aia3a4  . 

Since  the  gain  factors  are  nonzero,  this  determinant  can  never  equal 
zero.  It  can  be  shown  that  for  a  generalized  network,  any  selection 
of  columns  that  describes  a  tree  can  be  arranged  to  form  a  matrix 
with  diagonal  elements  equal  to  the  gain  factors  and  the  lower 
off-diagonal  elements  equal  to  zero.  Thus  | B | ,  the  determinant,  can 
always  be  written  as 


|  b|  -  +  (  n  a.  ) 

k£MB 
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where  M  is  the  set  of  arcs  in  the  basis  tree. 

D 

In  order  to  form  a  tree  in  the  manner  just  shown,  one  of  the 
arcs  chosen  must  originate  at  the  slack  node.  Otherwise  it  is 
impossible  to  arrange  the  arcs  so  that  the  first  column  has  only  one 
nonzero  element.  Pictorially,  we  will  represent  such  a  basis  as  a 
directed  tree  rooted  at  the  slack  node.  Where  it  is  necessary,  we 
will  use  mirror  arcs  to  obtain  the  directed  tree.  Figure  4.2 
illustrates  two  basis  trees  for  the  example  problem. 

Next,  consider  the  case  in  which  the  arcs  chosen  form  a 
cycle.  For  example,  select  arcs  3,  4,  5  from  Figure  4.1,  so  that 


and 


B 


-  a3a5. 


If  B  f  0  then  the  cycle  is  an  acceptable  basis.  We  note  that  the 
cycle  forms  a  basis  if  and  only  if 


This  has  a  graphical  interpretation  which  is  illustrated  in  Figure 
4.2.  Starting  at  node  2  and  passing  around  the  cycle,  we  encounter 
arcs  3,  5,  -4.  The  gain  of  the  cycle  is  then  defined  as  the  product 
of  the  gains  that  form  a  cycle.  We  assign  the  symbol  B  to  the  cycle 
gain.  Thus  for  the  example 

g  -  a3a5  /  a4 

and  this  cycle  would  be  an  acceptable  basis  if  B  does  not  equal  one. 

In  general,  it  can  be  shown  that  for  a  basis  component  with 
that  includes  cycle  arcs  Mg  (M 
the  associated  matrix  will  be  proportional  to 

11  \  -  n  a.  . 

keMB  ke  Mjj-Mg 

We  represent  a  basis  component  that  Includes  a  cycle  as  a  directed 
graph  with  all  noncycle  arcs  directed  as  if  the  cycle  were  the  root. 
It  is  convenient  to  think  of  these  arcs  as  a  tree  rooted  at  a  cycle. 
A  particular  basis  may  be  a  combination  of  components  as  illustrated 


arcs  Mg 


c  C  Mg) *  the  determinant  of 


in  Figure  4.4 


4.7  Pointer  Representation  of  the  Basis 


To  represent  the  components  of  a  tree,  three  labels  are 
assigned  to  each  node.  For  node  i,  the  three  labels  would  be  the 
back  pointer  Pg(i),  the  forward  pointer  Pp(i),  and  the  right  pointer 
PR(i).  The  back  pointer  is  the  unique  arc  which  terminates  at  node 
i.  The  forward  pointer  is  the  terminal  node  of  the  arc  originating 
at  node  i.  The  right  pointer  is  the  node  appearing  to  the  right  of 
node  i  in  the  tree.  The  basis  can  also  be  represented  with  the 
preorder  traversal  method  (37).  In  this  method  the  preorder 
traversal  list  Pp  simply  stores  the  pointer  representation  where 
Pp(i)  is  set  equal  to  the  node  number  to  which  node  i  points.  In 
our  descriptions,  as  well  as  in  Jensen  and  Barnes  (46),  the  root 
node  is  always  the  slack  node.  Figure  4.5  illustrates  these 
concepts. 

4.8  Linear  Programming  Results  for  Network  Models 

In  order  to  relate  material  presented  in  Chapter  3  and  to 
help  lay  the  foundation  for  much  of  the  material  presented  in 
Chapters  5  and  6,  it  seems  appropriate  to  summarize  some  of  the 


major  results  of  linear  programming  theory  that  pertain  to  network 
flow  problems.  In  terms  of  the  parameters  and  terms  previously 
defined  In  this  chapter,  the  generalized  single  commodity  linear 
minimum  cost  flow  problem  may  be  presented  as  the  linear  programming 
primal  problem: 


m 

Min  E 

Vk 

(4.1) 

k-1 

S.T.  I 

fk  “  2  \fk  ”  bi 

i*l f • • • p 

(4.2) 

ke 

fk  <  ck 

Mqi  ke  “ti 

k“l  |  •  •  •  )  TD 

(4.3) 

fk>0 

k*  1  p  •  •  •  p  ni 

(4.4) 

where  node  n  Is  the  slack  node. 

The  dual  to  the  linear  program,  may  be  written  as 


n-1  m 

Min  E  tr^b^  ^  £  6  ^c^ 

(4.5) 

1-1  k-1 

S.T.  w1  -  ak1Tj  +  6k  >  ~hk 

k-1 ,. . . ,m 

i-O(k) ,  j-T(k) 

(4.6) 

■rr ^  unrestricted 


1-1 


n-1 


(4.7) 


The  primal-dual  conditions  for  an  optimal  solution  of  the 
linear  program  have  been  specialized  to  the  generalized  minimum  cost 
network  flow  problem  in  the  following  three  theorems. 

THEOREM  1 :  Given  a  solution  F  to  the  primal  problem  and  a  solution 
[  it  , 6  ]  to  the  dual  problem,  the  solutions  are  optimal  if  and  only 
if: 

1.  F  is  feasible  for  the  problem;  that  is 
(4. 2-4. 4)  are  satisfied. 

2.  [  tt  ,6  ]  is  feasible  for  the  dual  problem; 
that  is  (4.6)  and  (4.8)  are  satisfied. 

3.  Complementary  slackness  is  satisfied;  that 
is 

(a)  If  fk(i,j)  >  0  then  n  i  -  a^  ^  +6  k  « 

-v 

(b)  If  ffe  <  ck,  then  5  k  -  0. 

(c)  If  -  ak  it j  +6k  >  -hk,  then  ffc(i.j) 

*  0. 

(d)  If  <$k  >  0,  then  fk  -  ck> 

Since  C  >_  0,  a  restricted  condition  of  dual  feasibility  may  be 
succinctly  stated  as  follows: 


THEOREM  2:  If  [  tt  ,  6  J  is  an  optimal  solution  to  the  dual  problem, 
then  6  k  «  max  [  0,-h^  -  +  a^  n  ^  J. 

Theorem  2  allows  the  optimality  conditions  for  the  problem  to  be 
written  as  follows: 

THEOREM  3:  Given  a  solution  F  to  the  primal  problem  and  a  partial 
solution  ir  to  the  dual  problem,  the  solutions  are  optimal  to  their 
respective  problems  if  and  only  if  the  following  considerations  are 
satisfied: 

1.  Primal  feasibility. 

2.  6  k  -  max  [  0,-hj^  -  77  *  +  j  ] 

(restricted  dual  feasibility). 

3.  Complementary  slackness: 

(a)  Tr  i  "  \  17  j  “  f°r  0  <  fk  <  cfc. 

(b)  fk  ■  0  for  tt^  -  ak  tt  j  >  -1^. 

(c)  fk  ■  ck  for  it  ^  -  ak  iTj  <  -hk. 

These  results  lay  the  foundation  for  various  solution  procedures 
used  to  solve  generalized  network  flow  problems,  which  will  now  be 


presented 


4.9  Solution  Algorithms  for  Network  Flow  Problems 

A  wide  variety  of  solution  algorithms  have  been  introduced 
to  solve  network  problems  of  various  types.  In  general,  they  are 
finite  iterative  procedures  designed  to  obtain  a  solution  that 
satisfies  the  conditions  stated  in  section  4.8:  primal  feasibility, 
restricted  dual  feasibility,  and  complementary  slackness.  The  major 
difference  between  the  algorithms  is  the  order  in  which  the 
conditions  are  satisfied.  Each  of  the  algorithms  generally  takes  on 
a  different  form  as  it  is  applied  to  various  problem  classes,  but 
the  steps  noted  below  are  consistently  followed.  While  the 
conditions  for  optimality  do  not  require  that  the  solutions  be 
basic,  maintaining  a  basis  will  often  result  in  a  savings  of 
computational  t ime . 

(a)  Primal  Approach  -  The  primal  approach 

iteratively  derives  F  and  it  such  that  F  is 
primal  feasible  while  attempting  to  achieve 
complementary  slackness.  If  we  define  e^, 
a  measure  of  the  violation  of  complementary 
slackness  for  arc  k,  as 


eck  -  max  l  dkfk’  <fk  •  ck>dk  1 

for  feasible  flow,  e^  will  be  positive  if 
and  only  if  complementary  slackness  is 
violated.  Let 


m 

Ec  "  Z  eck* 

C  k-1  C 

When  Ec  is  driven  to  zero,  F  is  optimal. 

ALGORITHM 

1.  Find  F  and  it  that  satisfy  primal 
feasibility. 

2.  Find  an  arc  such  that  e^  >0.  If 
there  are  none,  stop.  Otherwise  go  to 
step  3. 

3.  Find  a  new  F  and  w  that  reduce  e  .  for 

ck 

the  arc  while  maintaining  primal 
feasibility.  Go  to  step  2. 


(b)  Dual  Node  Infeasible  Algorithm  -  The  dual 
node  infeasible  algorithm  maintains 


complementary  slackness  and  satisfies  all 

primal  feasibility  requirements  except 

/ 

conservation  of  flow.  Let  b^  be  the 
external  flow  requirement  at  node  i  that 
would  satisfy  conservation  of  internal 
flows  at  node  i,  under  current  arc  flow 
requirements.  Let  e^  ■  |  b^  -  b^  |  be  a 
measure  of  the  infeasibility  for  node  i  and 
let 


N-l 

EN  *  eNi  * 

In  this  approach  we  achieve  optimality  by 
iteratively  changing  F  and  tt  in  such  a  way 
to  ultimately  force  to  zero. 


ALGORITHM 

1.  Find  ^  and  F  £  0  that  satisfy 

complementary  slackness  and  the  arc 

/ 

capacity  restrictions.  Let  b.  be  the 


node  flows  required  to  obtain 
conservation  of  flow. 

2.  Find  a  node  such  that  e,T.  >0.  If 

N1  — 

there  are  none,  stop.  Otherwise  go  to 
step  3. 

3.  Find  a  new  F  and  it  that  reduces  the 
infeasibility  of  the  node  while 
maintaining  complementary  slackness 
and  arc  feasibility.  Go  to  step  2. 

c)  Dual  Arc  Infeasible  Algorithm  -  During  each 
Iteration  of  the  dual  arc  infeasible 
approach,  complementary  slackness  and 
conservation  of  flow  conditions  are 
satisfied.  However,  one  or  more  arcs  will 
have  flows  above  their  capacities  or  below 
zero.  Let  e^  -  max  [  ~ffc,  0,fk  -  ck  ]  be 
a  measure  of  infeasibility  of  arc  k.  Let 


measure  the  infeasibili'ty  of  the  network. 

The  dual  arc  Infeasible  algorithm  achieves 

optimality  by  iteratively  modifying  F  and  tt 

to  force  E,  *  0. 

A 

ALGORITHM 

1.  Find  an  initial  F  and  tt  that  satisfies 
complementary  slackness  and 
conservation  of  flow. 

2.  Find  an  arc  k  for  which  e^  >0.  If 
none  exist,  stop.  Otherwise  go  to 
step  3. 

3.  Modify  F  and  ir  in  order  to  reduce  e^ 
while  maintaining  conservation  of  flow 
and  complementary  slackness.  Go  to 
step  2. 

4. 10  Network  Manipulation  Subroutines 

At  this  point  we  have  presented  all  the  necessary  background 
except  for  one  major  area.  In  Chapters  3  and  6,  numerous  references 
will  be  made  to  various  subroutines  developed  by  Jensen  and  Barnes 


(46).  It  would  be  impractical  to  recreate  the  development  of  these 


subroutines  at  this  time.  Instead  the  following  material  is 
review  of  the  purpose  of  each  subroutine  referenced. 


a)  ADDTRE-  To  add  an  arc  k(i,j)  to  a  forest  of 

trees.  Node  i  and  j  must  be  in  different 
trees  and  node  j  must  be  the  root  of  a 
tree. 

b)  CYCLE-  To  compute  the  gain  of  a  cycle  defined 

by  the  backpointers  and  the  cost  of 
circulating  one  unit  of  flow  around  the 
cycle. 

c)  DELTRE-  To  delete  an  arc  from  the  pointer 

representation  of  the  tree. 

d)  DSHRTG-  To  derive  the  shortest  path  from  node 

s  to  node  t  in  a  generalized  network  when 
all  costs  are  positive  and  all  gains  are  £ 
1. 

e)  DUAL-  To  compute  the  dual  variables  for  the 

basis  network  rooted  at  a  given  node  II. 

f)  FLOWG-  To  change  the  flow  in  each  arc  by  an 

amount  prescribed  by  MF  and  G(j).  (See 
Chapter  5  for  a  definition  of  these 
quantities) . 

g)  ORIG-  To  determine  the  list  of  arcs 

originating  at  node  I. 

h)  ORIGSG-  To  accept  an  arc  data  item  and  store 

it  in  an  arc  list  ordered  by  ascending 
origin  node. 

i)  PSHRTG-  To  solve  the  shortest  path  problem 

for  the  generalized  network. 

j)  READG-  To  read  and  store  node  and  arc  data 

for  the  generalized  minimum  cost  flow 
problem. 


k)  ROOTG-  To  find  the  list  of  arcs  (LISA)  and 

the  list  of  nodes  (LISN)  that  are  in  the 
directed  tree  rooted  at  node  IROOT. 

l)  TERM-  To  determine  the  list  of  arcs 

terminating  at  node  I. 

m)  TERMS-  To  create  LT,  a  list  of  arc  indices  in 

order  of  increasing  terminal  node. 

n)  TREINT-  To  construct  a  pointer  representation 

of  a  tree  given  knowledge  of  the  back 
pointer  arcs  which  make  up  the  tree. 

o)  TRECHG-  To  delete  an  arc  k  from  the  basis 

tree,  insert  another  arc  k  in  the  basis 
tree,  and  redirect  certain  arcs  in  the  tree 
to  maintain  a  directed  tree. 


CHAPTER  5 


PARAMETRIC  ANALYSIS  OF  A  GENERALIZED  NETWORK 


5.1  Introduction 


In  this  chapter  we  will  describe  an  algorithm  to  perform 
parametric  analysis  on  the  capacity  vector  and  the  external  flow 
vector  of  a  generalized  network.  The  mathematical  development  of 
this  material  as  related  to  linear  programming  theory  can  be  found 
in  Chapter  3,  sections  3.8  and  3.9.  We  will  first  develop  the 
detailed  algorithms  for  the  variation  of  the  capacity  vector  and 
then  extend  this  to  the  external  flow  vector. 

5.2  Parametric  Analysis  for  Arc  Capacities 

The  parametric  sensitivity  analysis  algorithm  described  in 
this  section  is  applicable  to  the  generalized  minimum  cost  flow 
problem  in  which  the  signs  of  the  arc  costs  and  arc  gains  are 
unrestricted.  The  algorithm  is  a  dual  algorithm  since  the  iterative 
step  first  selects  an  arc  to  be  deleted  from  the  basis  and  then 
selects  another  arc  to  be  added  to  the  basis.  We  may  state  the 
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general  steps  of  the  algorithm  as  follows: 


1*  Solve  the  original  problem  to  optimality. 

2.  Specify  the  list  of  candidate  arcs  and  the 
list  of  parametric  parameters. 

3.  Find  the  flow  augmenting  trail  or  trails 
which  are  generated  by  nonbasic  candidate 
arcs  .  Determine  the  arc  to  leave  the 
basis. 

4.  Find  an  arc  to  enter  the  basis. 

5.  Change  the  basis  by  deleting  the  leaving 
arc  and  inserting  the  entering  arc.  Modify 
the  dual  variables  to  satisfy  complementary 
slackness  for  the  new  basis.  Return  to 
step  3. 


Before  stating  the  algorithm  in  detail,  each  step  of  the  outline  is 
discussed. 


5.3  Initial  Optimal  Solution 

The  initial  problem  can  be  solved  to  optimality  by  either  a 
primal  simplex  algorithm  or  a  dual  flow  augmentation  algorithm  for 
the  generalized  minimum  cost  flow  problem.  For  details  of  these 
solution  methods,  which  will  not  be  repeated  here,  see  PGAINS  [46, 
pg  332]  and  INCREMG  [46,  pg  314] . 


(c^)  are  shown  as  dashed  lines.  The  following  sections  describe  how 
the  marginal  flow  change  in  each  basic  arc  is  determined  for  various 
changes  in  the  network  parameters.  Throughout  this  discussion,  the 
networks  of  Figures  5.1,  5.2,  and  5.3  are  used  as  examples. 

5.5.1  Flow  Change  at  a  Node 

If  the  external  flow  changes  at  a  node,  the  marginal  changes 
in  the  basic  flows  can  be  predicted  by  the  procedures  of  this 
section.  The  same  procedures  will  be  expanded  in  the  next  section 
to  find  the  marginal  flow  changes  caused  by  capacity  changes  of  a 
nonbaslc  arc. 

In  order  to  develop  the  relationships  that  guide  the  flow 
changes  at  a  node,  we  require  some  additional  definitions.  A  flow 
augmenting  trail  is  a  list  of  basic  arcs.  The  first  arc  on  the 
trail  may  either  originate  at  the  slack  node,  node  n,  or  at  a  node 
on  a  cycle.  The  last  arc  on  the  trail  may  terminate  at  any 
arbrltary  node  v.  Later  we  will  specify  exactly  what  this  last  node 
will  be  in  our  applications.  We  will  define  the  node  gain  for  any 
node  u,  Yu,  u£  N,  as  the  inverse  product  of  the  arc  gains  on  the 
trail  from  node  u  to  the  last  node  on  the  trail.  We  will  denote  T^ 
as  the  set  of  arcs  on  the  trail  from  node  u  to  the  last  node  on  the 


trail  and  Y 


1  if  u  is  the  last  node  on  the  trail 


Once  the 


value  of  Yu  is  computed,  where  u  is  any  node  in  the  network,  we  can 
compute  the  flow  in  a  basic  arc  k(i,j)  as  a  function  of  the  flow 
change  at  node  J .  Let  A  be  the  flow  increase  at  node  j .  The 
corresponding  increase  in  arc  k  is 

Ak  ’  YjA/ak. 

To  illustrate  the  concept  of  flow  augmenting  trails, 
consider  a  trail  from  node  11  to  node  9  in  Figure  5.1.  In  this  case 
the  arcs  on  the  trail  are  [14,6,2,22]  and  the  nodes  on  the  trail  are 
[9,6,3,1,11].  Note  the  trail  is  found  by  tracing  the  backpointers 
from  the  end  of  the  trail.  We  may  calculate  the  node  gains  for  each 
node  on  the  trail  as 

Y9  -  i 

Y6  *  1/a14 

Y  3  *  l/(ai4a6^ 

Y  j  -  l/(aj^a6a2) . 

y  11  "  1/(a14a6a2a22)* 

As  an  example  of  a  trail  including  a  cycle,  consider  the  trail 
starting  at  node  6  and  goes  to  node  7  in  Figure  5.2.  The  arcs  on 
the  trail  are  [11,5,7,-6,12]  and  the  nodes  are  [7, 5, 2, 3, 6].  Again 


we  note  that  the  trail  found  by  tracing  backwards  from  node  7  by 
using  backpointers.  We  may  calculate  the  node  gains  for  each  of  the 
nodes  on  the  trail  as 


Y5  *  U/anK  3  /  6  -  1) 

Y2  -  [l/(ana5)](B  /  6  -  1) 

Y3  «  [l/(a11a5a7)](6  /B  -  1) 

Yfe  »  [l/(a11a5a7a_6)]( 3  /3  —  1). 

where  g  is  the  cycle  gain  described  in  section  4.6. 

Thus  in  general,  If  node  u  is  not  on  a  cycle  we  may  write 


Y 

u 


U.  11  _  *k 
k£  T 


and  if  node  u  is  on  a  cycle 


(5.1) 


Y  -  (1/  n  ak)(  S/6  -  1) 
ke  T 

u 


and  Y  »  1  if  u  is  the  last  node  on  the  trail. 


(5.2) 


The  effect  of  a  capacity  change  in  a  nonbasic  arc  depends  on 

whether  the  flow  on  this  arc  is  zero  or  at  capacity.  If  the  flow  is 

zero,  there  is  no  effect  and  no  analysis  is  necessary.  If  the  flow 

is  at  capacity,  then  as  the  capacity  changes  the  flow  in  the  arc 

must  change  simultaneously  so  that  the  arc  flow  remains  at  capacity. 

The  arc  will  remain  nonbasic  with  flow  at  its  upper  bound.  The 

basic  flows  must  change  to  accommodate  the  flow  change  in  the 
« 

nonbasic  arc.  Increasing  or  decreasing  the  flow  in  the  nonbasic  arc 
x^  has  the  effect  of  changing  the  node  flows  at  its  two  terminal 
nodes. 

For  a  nonbasic  arc  x^  two  augmenting  trails  are  defined. 
There  is  an  augemtlng  trail  that  starts  at  the  slack  node  n  or  at  a 
cycle  and  ends  at  the  origin  node  of  x^.  The  second  flow 
augmenting  trail  starts  at  the  slack  node  n  or  at  a  cycle  and  ends 
at  the  terminal  node  of  x^.  These  two  flow  augmenting  trails  may 
have  no  nodes  in  common,  or  one  or  more  •"'des  in  common.  There  are 
many  possible  combinations  of  flow  augmenting  trails.  Consider  the 
nonbasic  arc  in  Figure  5.3(a).  In  this  case  both  flow  augmenting 
trails  originate  at  a  cycle  and  end  at  the  origin  node  and  terminal 
node  of  x^  respectively.  We  should  note  that  the  trails  have  no 
nodes  in  common.  In  Figure  5.1(b),  nonbasic  arc  20  forms  two 


trails.  One  trail  starts  at  node  11  and  ends  at  node  9.  It  has  arcs 

[14,6,2,22]  and  nodes  [9,6,3,1,11].  The  second  trail  begins  at  node 

11  and  ends  at  node  10.  It  has  arcs  [15,11,5,7,2,22]  and  nodes 

[10,7,5,2,3,1,11].  In  this  case  the  trails  intersect  at  node  3  and 

have  common  arcs  2  and  22  and  common  nodes  3,  1,  and  11. 

In  Figure  5.2  nonbasic  arc  15  also  forms  two  trails.  The 

first  trail  begins  at  node  5  and  ends  at  node  7.  The  trail  contains 

arcs  [11,5,7,-6,12]  and  nodes  (7, 5, 2, 3, 6].  In  this  case  the  trail 

originates  at  a  cycle.  The  second  trail  begins  at  node  6  and  ends 

at  node  10.  It  contains  arcs  [20, -19, -9, 12,5,7 , -6]  and  nodes 

[10,9,4,6,5,2,3].  In  this  example  both  trails  have  a  common  cycle. 

In  section  5.5.2  we  defined  Ty  as  the  set  of  arcs  on  the 

trail  from  u  to  the  last  node  on  the  trail  and  Y  *  1  if  u  is  the 

u 

last  node  on  the  trail.  For  a  nonbasic  arc  x^d.j)  for  which  the 
capacity  is  to  be  varied,  this  is  no  longer  true.  We  must  adjust 
the  Y  's  at  the  origin  and  the  terminal  nodes  of  x^d.j)  to  account 
for  the  effects  of  the  parametric  change  we  are  going  to  make  on  arc 
x^d.j).  In  particular  for  the  origin  node  of  arc  x^i.j),  node  i, 
which  is  the  last  node  of  one  trail,  Y  a  p^.  For  the  terminal 
node  of  arc  x^,  node  j,  which  is  the  last  node  on  the  second  trail, 

yi  m  -PkV 

Each  trail  generated  by  a  nonbasic  arc  x^  will  create  node 
gains  for  each  node  on  the  trail.  The  following  examples  will 
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demonstrate  that  by  using  a  series  of  additive  operations,  Y  can 
be  computed  depending  on  where  node  u  is  located  on  the  trail. 

Example  1 :  Consider  the  network  shown  in  Figure  5.1(b).  As  noted 

earlier,  candidate  arc  20  forms  two  flow  augmenting  trails. 
Calculating  the  node  gains  for  the  flow  augmenting  trail  starting  at 
node  11  and  ending  at  node  9,  we  have 

Y9  "  p20 
Y6  ‘  p20/a14 

Y  3  *  p20/(a14a6) 

Y1  "  p20/(a14a6a2) 

Yll“  p20/(a14a6a2a22) ' 

However  we  must  also  consider  the  second  flow  augmenting  trail  from 
node  11  to  node  10.  In  this  case  we  calculate  the  node  gains  as 

Y10  "  ~p20a20 

Y  7  *  "p20a20/a15 

Y  5  “  _p20a20/(a15all) 

Y2  "  _p20a20/(a15aUa5) 

Y  3  “  _p20a20/(a15alla5a7) 
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Y1  *  "p20a20/(a15alla5a7a2) 

Y11  "  ■p20a20/(a15alla5a7a2a22)* 

We  note  that  nodes  3,  1,  and  11  have  been  encountered  twice,  once 
with  the  flow  augmenting  trail  ending  at  node  9  and  once  with  the 
flow  augmenting  trail  ending  at  node  10.  Thus  the  node  gains  for 
nodes  3,  1,  and  11  must  be  adjusted  for  these  nodes  on  a  common 
augmenting  trail.  Therefore 

Y  3  ”  p20/(a14a6)  "  p20a20/(a15alla5a7) 

Y1  "  p20/(a14a6a2)  '  p20a20/(a15alla5a7a2) 

Yll“  p20/(a14a6a2a22^  ~  p20a20^a15alla5a2a22)# 

Since  the  network  is  linear,  we  account  for  the  node  gains  on  the 
common  part  of  the  flow  augmenting  trail  by  adding  the  node  gain 
values  of  the  first  flow  augmenting  trail  to  the  node  gains  of  the 
second  flow  augmenting  trail.  The  node  gains  on  the  unique  part  of 
the  flow  augmenting  trails  remain  unchanged. 

Example  2;  Consider  the  network  shown  in  Figure  5.2.  We  will  now 
calculate  the  node  gains  for  the  flow  augmenting  trails  generated  by 
arc  15.  To  calculate  the  node  gains  for  the  flow  augmenting  trail 
starting  at  node  6  and  terminating  at  node  7  we  have  : 


AD-A145  558  PARAMETRIC  ANALYSIS  FOR  GENERALIZED  NETWORK  FLOW  2/m 

PROBLEMSOJ)  AIR  FORCE  INST  OF  TECH  HRIGHT-PATTERSON  AFB  * 
OH  M  E  BAUM  MAY  84  AFIT/CI/NR-84-38D 


UNCLASSIFIED 


F/G  12/1 


NL 


Y5  "  lp15/all]  (S  -  L) 
y  2  *  tp15/(alla5^  <S  /B  -  O 
Y3“  ^p15/(alla5a7^  <0/6  “  1) 

Y6”  ^pl5/(alla5a7a-6^  <B/B“  1). 

To  calculate  the  node  gains  for  the  flow  augmenting  trail  starting 
at  node  3  and  terminating  at  node  10 ,  we  have  : 


Y  10  *  “PlS^S 

Y  9  "  -p15a15/a20 

y4  *  -p15a15/(a20a-19) 

Y  6  [_p15a15/(a20a-19a-9)I  (S/8"  l> 

Y 5  *  [_p15a15/(a20a-19a-9a12)]  (6/8_  X> 

Y  2  “  ^~p15a15/(a20a-19a-9a12a5)1  (  S/B-  1) 

Y3  ■  I ~ P 1 5*15/ ^ a20a— 19a— 9a 1 2a5a7^ ^  ^  ® I  ® • 


We  note  that  nodes  2,  3,  5,  6  appear  twice.  Thus  these  node  gains 
must  be  adjusted  since  they  apppear  on  a  common  augmenting  trail. 
For  example  the  adjusted  node  gain  of  node  2  is 


and  Che  node  gains  for  nodes  3,  5,  and  6  may  be  adjusted  similarly 
The  node  gains  on  Che  unique  pares  of  Che  flow  augmenting  trails 
remain  unchanged. 


5.5.3  Capacity  Change  on  a  Set  of  Nonbasic  Arcs 

When  we  consider  more  than  one  nonbasic  arc,  we  must  take 
into  account  the  combined  effects  of  numerous  intersecting  flow 
augmenting  trails.  We  should  note  that  if  any  of  the  nonbasic  arcs 
have  zero  flow,  no  further  analysis  is  needed  for  this  arc.  The  arc 
will  remain  nonbasic  as  the  change  on  the  capacity  will  have  no 
effect  on  the  flow  since  the  flow  is  at  zero. 

If  we  consider  the  example  network  in  Figrure  5.2,  we  see 
that  there  are  two  nonbasic  arcs,  arcs  14  and  15.  For  each 
candidate  arc,  two  flow  augmenting  trails  are  produced.  For  arc  15, 
one  trail,  starts  at  a  cycle  and  ends  at  node  7  with  arcs 
[11,5,7,-6,12]  and  nodes  [7, 5, 2, 3, 6].  The  second  trail  created  by 
arc  15  has  arcs  [20,-19,-9,12,5,7,-61  and  nodes  [10,9,4,6,5,2,3]. 

For  arc  14  the  first  flow  augmenting  trail  has  arcs  [12, 5, 7, -6]  and 
nodes  [6, 5, 2, 3].  The  second  trail  created  by  arc  14  has  arcs 
f-19,-9,12,5,7,-6]  and  nodes  [9, 4, 6, 5, 2, 3] . 
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In  each  case  to  calculate  the  Y  values  for  each  node  on  -the 
trails,  we  see  that  certain  nodes  are  on  unique  trails,  while  other 
nodes  are  part  of  intersecting  trails.  For  example,  node  7  has  a 
unique  Y  value  as  no  other  flow  augmenting  trails  pass  through  this 
node.  In  particular 

y7  "  PI5* 

This  is  contrasted  to  node  6  which  is  on  flow  augmenting 
trails  which  end  at  nodes  7,  10,  and  9  respectively.  Thus  Y  must 
take  into  account  the  additive  contributions  from  each  of  the  trails 
which  pass  through  node  6.  In  particular 

Y6  “  ^p14  "  p14a14^*-19a-9^  +  p15/(alla5a7a-6) 
"p15a15/(a20a-19a-9^ 

In  order  to  write  the  general  formulas  to  calculate  the  node  gains, 
we  need  some  additional  notation.  Define  Tu  to  be  the  set  of  flow 
augmenting  trails  that  pass  through  node  u.  Let  T  be  the  set 

U,T 

of  arcs  on  the  trail  TC^U  that  originate  at  node  u  and  go  to  the 
end  of  the  trail.  Then  we  can  write  the  general  formulas  for  the 
calculation  of  the  node  gains  as 


is  Increased,  there  will  be  no  change  in  any  of  the  flows  in  the 
basic  arcs  as  the  candidate  arc  under  consideration  will  remain 
nonbasic  with  zero  flow.  If  the  capacity  is  being  reduced  no  effect 
will  occur.  However,  we  should  note  that  when  the  candidate  arc's 
capacity  becomes  equal  to  zero,  the  candidate  arc  is  nonbasic  at  its 
upper  bound  and  that  any  further  reductions  in  the  candidate  arc's 
capacity  will  cause  the  problem  to  become  Infeasible. 

Lastly,  we  must  consider  nonbasic  candidate  arcs  which  have 
flow  equal  to  their  upper  bound.  Changes  in  the  capacity  by  either 
increasing  or  decreasing  the  nonbasic  arc  under  consideration  will 
cause  flow  in  the  arc  to  change  simultaneously.  These  arcs  will 
remain  nonbasic  at  their  upper  bounds.  However,  we  must  adjust  the 
basic  flows  to  accommodate  these  flow  changes  in  the  nonbasic 
candidate  arcs  at  their  upper  bound  in  the  manner  of  section  5.5.3. 

5.5.5  Computer  Implementation  of  Flow  Augmenting  Trails 

The  subroutines  used  to  calculate  the  node  gains  for  the 
flow  augmenting  trails  are  PATH,  TRAIL2,  and  CYCLEG.  Subroutine 
PATH  is  called  by  the  main  program  PSENG  when  a  flow  augmenting 
trail  is  to  be  calculated  for  a  nonbasic  candidate  arc  with  flow 
equal  to  its  upper  bound.  TRAIL2  is  called  by  PATH  when  a  node  is 
encountered  on  a  common  part  of  a  flow  augementing  trail.  CYCLEG  is 


used  to  compute  the  gain  of  any  of  the  cycle(s)  found  during  PATH 
and  to  update  the  node  gain  for  each  of  the  nodes  on  the  cycle(s). 

For  Che  computer  implementation  of  the  concepts  developed  in 
section  5.5.1  through  5.5.4,  we  need  to  define  the  variables  which 
are  used  in  the  BASIC  code.  For  each  node  on  the  flow  augmenting 
trail(s)  we  assign  a  node  check  value,  denoted  IK(IJ),  to  enable  us 
to  determine  if  a  node  has  been  previously  encountered  on  a  flow 
augmenting  trail.  We  use  the  computer  coding  technique  to  set  the 
value  of  MULT  to  equal  +1  if  we  are  begining  the  trace  of  a  flow 
augmenting  trail  that  terminates  at  the  origin  node  of  x^,  or  we  set 
MULT  equal  to  -1,  if  we  are  starting  a  trace  of  a  flow  augmenting 
trail  that  ends  at  the  terminal  node  of  x^.  LN  refers  to  the  list 
of  nodes  encountered  and  LA  refers  to  the  list  of  arcs  encountered 
on  the  trail(s).  IC  is  a  counter  to  indicate  the  number  of  nodes  we 
have  encountered  and  CK  is  an  indicator  if  this  node  has  only  been 
encountered  on  this  trail  (ie,  on  a  cycle).  The  RC(Z)  list  is  a 
list  of  root  nodes  of  all  the  cycles  encountered  during  the  various 
passes  through  the  flow  augmenting  trall(s).  The  G(IJ)  list 
contains  the  node  gain  values  for  each  node  encountered  on  the  flow 
augmenting  trail. 

We  should  note  that  the  subroutine  PATH  in  conjunction  with 
TRAIL2  and  CYCLEG  differ  from  the  subroutine  PATHPG  [46,  pg.  3231  in 
three  ways.  First,  for  each  candidate  arc  x^,  the  parametric 
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parameter  must  be  taken  into  consideration.  Second,  the  Y  ' a  are 
calculated  additively  for  each  flow  augmenting  trail  created  by  a 
ca  didate  arc.  Finally,  if  a  cycle  is  encountered,  the  root  node  of 
the  cycle  is  entered  on  the  root  cycle  list  [RC(Z>].  Then  after  all 
the  flow  augmenting  trails  have  been  generated,  the  cycle  gains  are 
calculated  for  each  cycle  in  the  root  cycle  list  and  then  each  node 
gain  on  a  particular  cycle  is  updated. 

As  an  example  of  the  use  of  PATH,  consider  the  application 
to  the  network  shown  in  Figure  5.2.  The  following  information  is 
known  concerning  the  arcs  whose  capacities  are  being  varied.  The 
candidate  arcs  are 

X  -  [  5,6,7,11,14,15,20] 

and  the  parametric  values  are 


P  -  [-. 3,1. 2,-2. 1,-0. 67, 0.1, -3,1. 4] . 


We  note  that  only  arcs  14  and  15  are  considered  in  this  analysis  as 
they  are  the  only  nonbasic  candidate  arcs  at  their  upper  bounds. 

The  list  of  basic  arcs  at  the  current  optimal  solution  is 
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PURPOSE:  To  find 
che  flow  augment¬ 
ing  trail  for  a 
generalized  net¬ 
work  and  compute 
T’s  for  each  node 
on  the  trail. 

1. (INTIAL)  If  node 
1J  has  never  been 
encountered  (IK( 
IJ)-0),put  IJ  in  LN 
and  calculate  7  for 
node  IJ.  Otherwise 
node  IJ  was  visited 
before.  Update  node 
gains  with  TRAIL2 
and  RETURN. 

2.  (ROUTE)  Identify 
arcs  on  the  flow 
augmenting  trail. 
Find  the  back 
pointer  of  the  cur¬ 
rent  node.  If  zero 
the  slack  node  has 
been  visited.  If 
not  find  next  node 
in  trail  and  check 
if  visited.  If  not 
add  IJ  to  LN  and  K 
to  LA.  Update  node 
gain  and  node 
check  value  ,  go  to 
2.  Otherwise, check 
if  node  only  visit- 
on  this  path  (IK  > 
CK) .If  so  put  IJ  on 
root  cycle  list  and 
return.  Otherwise, 
IJ  has  been  visited 
on  a  previous  flow 
augmenting  trail. 
Call  TRAIL2  to  up¬ 
date  7's  on  common 
trail  and  return. 


CPATB(IJ,MULT,KE,XK,P,IK,CK)‘) 


INITIAL 

IK(IJ)  -  0  /l 

LN(IC+1):-IJ 

KE*MULT  >  0 

KE*MULT  >  0 

G:-P*MULT  G:-P*A(XK)‘* 
MULT 

G:-P*MULT  P:»P*A(XK)* 
MULT 

G(IJ):-  G(IJ)  +  G 

(  TRAIL2(IJ,G) ) 

- >ROUTE 

RETURN 

Subroutine  PATH 
Figure  5.4 
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PURPOSE:  Update 
values  on  cotsaon 
trail  found  In 

PATH.  INITIAL 


(  TRAIL2(IJ7gT~) 


1 .  (INTIAL)  Set  all 
node  check  values 
to  zero.  Update 

y  value  of  coaaon 
node  IJ. 

2. (ITBACK)  Look  at 
back  pointer  of  IJ. 
If  It  Is  zero,  the 
slack  node  has  been 
encountered,  re¬ 
turn.  If  not  check 
If  node  has  been 
visited  before,  if 
not  update  check 
value  and  7  for 
node.  Go  to  2.  If 
node  has  been  vis¬ 
ited  before,  a  root 
of  a  cycle  has  been 
found,  return. 


!I:-IJ,  CI:-0 

u 

U 

Q:-  1  TO  N 

IL(Q) :■  0 

[L(HI):*  1,  G(HI):--C(F 

[I)+G 

ITBACK 


K:-  PB(IJ) 

K  -  0 

K 

K 

>  0 

_ A 

IJ:-0(K) 

IJ:-T(-K) 

CI:-CI+l  | 

IL(IJ)  >  0  /N 

\ 

K  ) 

>  0  y/N 

G:-G*/A(K) 

G:«G*A(-K) 

G(IJ):-G(IJ)+G 

RETURN 

- MTBACK 

Subroutine  TRAIL2 
Figure  5.5 


PURPOSE:  Compute 
the  gain  of  a  cycle 
defined  by  the  back 
pointers  and  to  up¬ 
date  the  node  gains 
for  each  node  on  a 
cycle  by  its  cycle 
gain. 

1.  (INITIAL)  Set  II 
to  be  root  node  of 
cycle  from  cycle 
list.  If  all 
cycles  have  been 
updated  (I  >  Z), 
return.  Otherwise 
calculate  cycle 
gain  for  cycle 
RC(I)  by  calling 
CYCLE. 

2.  (MULTIPLIER)  Up¬ 
date  y  values  for 
for  each  node  on 
cycle  RC(I).  Up¬ 
date  I  and  go  to  1. 


Subroutine  CYCLEG 
Figure  5.6 


4 


The  basic  arcs  are  shown  in  Figure  5.2  as  solid  lines  and  the 
candidate  arcs  at  capacity  by  dashed-dot  lines. 

At  the  conclusion  of  PATH  for  the  candidate  arc  14  ending  at  node  6, 
the  following  information  has  been  obtained: 

Node  List:  LN  -  [6,5,2, 3} 

Arc  List  :  LA  -  [12, 5, 7, -6] 

Nodes  Visited:  IK  -  [0,3,4,0,2,1,0,0,0,0,01 

Node  Gains:  G  -  [0,0.057,0.067,0,0.143,0.1,0,0,0,0,0] 

Cycles  Encountered:  Z  -  1 

Root  Cycle  List  :  RC( 1)-  6 

Nodes  Encountered:  CK  -  4. 


We  have  simply  followed  the  back  pointers  from  node  1^  -  6 
until  a  node  is  encountered  that  has  been  visited  before.  In  this 
case  the  node  is  node  6,  since  a  cycle  is  formed.  The  list  IK 
serves  the  purpose  to  mark  nodes  that  have  been  visited.  At  this 
point  node  6  is  entered  on  the  root  cycle  list  [RC(1)  -  6]  and  CK  is 
sec  Co  4  to  indicate  a  cycle  has  been  found  on  this  augmenting  path. 
LN  contains  the  nodes  in  the  order  that  they  were  encountered  by 
PATH,  and  LA  contains  the  backpointers  of  the  nodes  in  LN.  These 
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backpointers  comprise  the  list  of  arcs  in  the  trail.  IC-is  a 

counter  to  indicate  how  many  nodes  have  been  encountered.  The  G 

list  indicates  the  value  of Y  for  each  node. 

u 

At  the  completion  of  PATH  for  arc  14  ending  at  node  9,  the 
lists  have  been  updated  to  read: 


LN  -  [6,5,2,3,9,41 
LA  -  [12, 5, 7, -6, -19, -9] 

IK  -  [0,3,4,6,2,1,0,0,5,0,01 
G  -  [0,-. 57, -.06, -.25, -.14, -.1,0,0, -.1,0, 01 
CK  «  6  Z  -  1  RC(1)  -  6. 

Note  that  when  node  6  is  encountered  IK(6)  -  1  and  CK  ■  6 
indicating  this  node  had  been  encountered  on  a  previous  cycle  and 
the  Yu's  need  only  be  updated  for  the  common  part  of  the  flow 
augmenting  path.  Recall  the  effects  of  the  cycle  gain,  6  /(8  -  l), 
will  be  accounted  for  after  the  calculations  for  all  the  flow 
augmenting  paths  are  completed. 

We  must  now  also  add  in  the  effects  for  the  candidate  arc, 

■  15.  After  the  completion  of  PATH  ending  at  node  10,  the  lists 
have  been  updated  again  to  be: 
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LN  -  [6,5,2,3,9,4,10] 

LA  -  [12,5,7,-6,-19,-9,20] 

IK  -  10,3,4,6,2,1,0,0,5,7,0] 

G  -  [0,2.85,3.36,6.12,7.14,5,0,0,2.45,2.55,0] 

CK  =■  7  Z  -  1  RC(  1 )  -  6. 

In  ROUTE  node  9  is  again  encountered,  so  the  additive  effects  for 
the  common  augmenting  path  updates  the  G  list. 

After  the  completion  of  PATH  ending  at  node  7,  the  lists 
have  been  updated  to  read: 

LN  -  [6,5,2,3,9,4,10,7] 

LA  -  [12,5,7,-6,-19,-9,20,11] 

IK  -  [0,3,4,6,2,1,8,0,5,7,01 

G  -  [0,1.65,1.94,6.12,4.14,3.72,-3,0,2.45,2.55,0] 

CK  -  9  Z  «  1  RC(1)  •  6. 

As  there  are  no  more  candidate  arcs  to  be  considered,  the 
cycle  gain  for  the  cycle  encountered  during  PATH  needs  to  be 
calculated.  In  this  case  it  is  rooted  at  node  6  and  contains  the 
nodes  2,  3,  5,  and  6.  The  cycle  gain  is  1.65  and  the  updated  G  list 


looks  like: 
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G  -  [0,4.19,4.93,6.12,10.49,9.44,-3,0,2.45,2.55,0] 


and  Che  lists  LN,  LA,  and  IK  remain  unchanged. 


5.6  Computing  Maximum  Parametric  Variation  of  Capacities  and 


Finding  the  Arc  to  Leave  the  Basis 


As  in  the  previous  sections,  we  assume  we  have  an  optimal 
solution  to  a  given  network  flow  problem.  We  define  the  maximum 
reference  number,  r,  as  the  maximum  flow  augmentation  that  will 
drive  a  basic  flow  to  zero  or  capacity  or  drive  a  candidate  arc's 
upper  bound  to  zero.  Then  as  we  change  the  capacity  on  a  candidate 


arc  x^,  the  new  capacity  on  the  candidate  arc  will  be 
ck  “  ck  +  Pkr 


(5.5) 


where  keX  and  c^  is  the  capacity,  is  the  parametric  parameter, 
and  r  is  the  reference  number  of  the  candidate  arc. 

We  now  need  to  determine  the  maximum  positive  increment  rfc 
that  will  drive  a  basic  flow  to  a  bound  or  cause  the  capacity  on  a 
candidate  arc  to  go  to  zero.  Assume  there  are  initial  flows  in  the 
arcs  on  the  flow  augmentation  trail(s),  identified  as  f^  or  f  k 
depending  on  whether  k  is  a  forward  or  a  mirror  arc.  We  should  note 


that  since  arc  k  is  basic  and  if  is  not  equal  to  zero,  the 


capacity  on  the  basic  arc  is  changing  as  well  as  the  flow.  Thus  we 
must  monitor  the  simultaneous  effects  of  both  changes.  Recalling 
that  is  defined  as  the  maximum  flow  change  that  will  drive  the 
flow  in  arc  k  to  one  of  its  bounds,  we  need  to  consider  two  cases. 

CASE  1:  k  is  a  forward  arc 

If  k  >  0,  then  the  flow  change  in  the  original  network 
implied  by  a  flow  increment  r  is 

Ak- Y/ak 

where  Y^  is  the  node  gain  at  the  terminal  node  of  arc  k.  Since  the 
initial  flow  on  arc  k  is  defined  by  f^,  then  the  flow  in  arc  k  after 
a  flow  increment  is 

fk*  4k-4k*  Y'v 

Since  r  is  greater  than  or  equal  to  zero,  the  direction  of  the  flow 
change  depends  only  on  the  sign  of  Y^/a^. 

SUBCASE  A  Yj/a^  ^0  :  WithYj/a^  positive,  the  arc  flow 
Increases.  Since  r  is  the  value  of  the  reference  number  that  will 
drive  a  basic  flow  to  zero  or  a  new  capacity,  c^  +  p^r,  the  new  flow 


will  be  limited  by 


fk  +  (  Yjr/ak>  <  ck  +  Pkr 

or  the  value  of  r  that  drives  the  basic  arc  k  to  its  upper  bound  is 
rk  *  (ck  "  fk)ak/(  Yj  "  PkV*  (5,6) 


We  should  note  that  if  -  p^a^  0,  the  new  capacity  is 
increasing  quicker  than  the  flow  is  changing.  This  means  this  arc 
would  remain  basic  since  the  new  flow  will  always  remain  less  than 
its  new  capacity.  We  set  r^  ■  00  ,  since  it  will  not  create  a 
limiting  value. 


SUBCASE  B  Yj/ak  <  0  •  In  this  case,  since  the  factor  Y^/a^  ls 
negative,  the  flow  in  arc  k  is  decreasing  with  r.  The  new  flow  will 
be  limited  by 


0  <  fk  .  <  ck  +  pkr. 

However,  two  possibilities  can  occur.  If  1  -  p^a^  <  0,  the  flow 

is  going  to  zero  quicker  than  the  capacity  and  therefore  the 
limiting  value  of  r  for  arc  k  is 


or  the  capacity  could  be  decreasing  faster  than  the  flow  is 


decreasing,  (  -  p^)  >0,  so  in  this  case 

rk  *  (ck  ‘  fk)ak/(Yj  -  PkV*  <5.8) 

CASE  2:  k  is  a  mirror  arc 

If  k  <  0,  then  the  flow  change  in  the  original  network 
Implied  by  a  flow  Increment  r  is 


where  j  -  o(-k).  Since  the  initial  flow  on  arc  k  is  defined  by  f  ,  , 

-k 

then  the  flow  on  arc  k  after  a  flow  change  is 


SUBCASE  A  ^j<0:  With^j<0,  the  arc  flow  will 
increase  and  the  value  of  r  which  drives  it  to  the  new  capacity, 
+  P_fcr»  will  be  limited  by 


f-k  “V  i  c-k +  p-kr 

or  the  limiting  value  when  Y  +  p_fc  <  0  is 
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If  Yj  +  p_k  0  the  new  flow  is  increasing  at  a  slower  rate  than 
the  capacity.  This  means  this  arc  would  remain  basic  since  the  new 
flow  will  be  less  than  the  new  capacity.  Therefore  we  set  rfe  »  00 
as  it  will  not  create  limiting  value. 

SUBCASE  B  Yj  >  0  :  In  this  case,  since  the  factor  Y^  >  0,  -Y^r 
is  negative,  and  the  flow  in  arc  k  is  decreasing.  Thus  the  flow 
will  be  limited  by 

0<f-k-  y  <  c* +  p-k'- 

Again  two  cases  exist.  If  Y ^  +  p_^  >  0,  the  flow  is  going  to  zero 
faster  than  the  capacity  and  thus 

rk  -  f_k/  Yj .  (5.10) 


Otherwise  the  capacity  could  be  decreasing  faster  than  the  flow  is 
decreasing,  so 


For  each  of  Che  six  possibilities  listed  above  in  equations 
(5.6)  through  (5.11),  r^  is  positive.  The  maximum  reference  number 
or  flow  increment  for  the  trail(s)  is  the  reference  number  that  will 
drive  the  flow  on  one  or  more  arcs  to  zero  or  its  capacity  without 
causing  any  other  flow  to  be  infeasible.  Thus 

Tj  -  Min  {rk}  (5.12) 

ke  Mq 

where  is  the  set  of  arcs  on  the  flow  augmenting  trail(s)  and  rfc 
is  determined  as  in  equations  (5.6)  through  (5.11). 

We  must  also  consider  the  possiblity  that  a  nonbaslc 
candidate  arc  will  have  its  capacity  to  zero.  Thus  for  each 
nonbaslc  candidate  arc  with  pk  <  0,  we  must  calculate 

rk  *  ck^~pk*  (5.13) 

If  a  nonbaslc  candidate  arc  has  a  positive  parametric  parameter,  p^ 

>  0,  we  set  r^  «  oo  a3  this  arc  will  not  create  a  limiting  value  of 
r. 

Again  we  need  to  calculate  the  maximum  reference  number  for 


the  nonbaslc  candidate  arcs 


Min  {r.  } 
k£XR 


(5.14) 


where  is  the  set  of  of  nonbasic  candidate  arcs  and  r^  is 
determined  as  in  equation  (5.13). 

We  must  pick  the  smallest  of  the  reference  numbers 
calculated  in  (5.12)  and  (5.14),  as  this  will  be  the  maximum 
reference  number  or  flow  increment  r  that  will  either  drive  a  basic 
arc  to  zero  or  capacity  or  a  nonbasic  arc's  capacity  to  zero.  Thus 


r  -  Min  {r^  r2>  . 


The  MFL0G2  subroutine  is  use  to  find  the  maximum  -flow  change 
in  the  trail,  given  the  Y  values  and  the  list  of  arcs  and  nodes  in 
the  trail(s).  PSENG  calls  MFL0G2  after  all  the  flow  augmenting 
trall(s)  have  been  scanned.  MFL0G2  will  determine  the  basic  arc 
that  limits  the  flow  change  on  the  flow  augmenting  trall(s)  and  thus 
the  arc  which  may  leave  the  basis.  As  in  Chapter  10  (46),  the  FLOWG 
subroutine  modifies  the  arc  flows  for  a  given  trail  or  trail(s)  and 
a  given  value  of  r. 

To  illustrate  some  of  the  types  of  calculations  that  are 
performed  in  MFL0G2 ,  consider  the  optimal  basis  shown  in  Figure  5.2. 
The  list  of  candidate  arcs  are 
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l 5, -6, 7 ,11,14,15,20] 

and  the  associated  parametric  parameters  are 


[-0.3, 1.2, -2. 1,-0. 67, 0.1, -3, 1.41. 


The  basic  arcs  on  the  flow  augmenting  trails  are 


[12,5,7,-6,-19,-9,20,11]. 


At  the  end  of  PATH  we  have  calculated  the  following  values  for 
then  nodes  1  through  11 


G  -  [0,4. 19, 4.93, ( 


As  an  example  of  the  calculations  preformed  in  MFL0G2  consider  arc 
20.  For  arc  20  we  have  \o^a20  *  ^.55  an<*  t*lis  Inean8  we  are 
case  1.  Calculating  r^Q,  we  have 

r20  “  ^c20“f20^a20^  Y10  '  p20a20^ 
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PURPOSE:  To  find 
che  maximum  flow 
augmentation  possi¬ 
ble  on  ell  basic 
arcs  In  LA  creat¬ 
ed  by  flow  augmen¬ 
ting  trails. 

1.  (INITIAL)  Set 
the  leaving  arc  to 
zero,  the  maximum 
flow  to  a  large 
number,  and  loop 
counter  to  zero. 

2.  (START)  Incre¬ 
ment  the  loop 
counter .  Check  if 
all  arcs  in  LA  have 
been  checked.  If 
so  return.  Other¬ 
wise  go  to  3. 

3.  (CALCULATE)  Set 
k  to  the  LA(I).  If 
k  Is  zero,  go  to  2. 
Otherwise  check  if 
k  Is  a  forward  arc. 
If  so  go  to  4, 
otherwise  go  to  3. 

4.  (FORWARD)  Set  j 
equal  to  the  term¬ 
inal  node  of  k.  De¬ 
termine  the  appro¬ 
priate  values  of 

D1  and  D2.  Go  to 

6. 


CALCULATE 


R:-LA(I),  Dl:-9999,  D2:-9999 

o 

■ 

u 

y/\ 

- >START 

K  >  0  yA 

- > FORWARD 

- >MIRROR 

FORWARD 


J:-T(K) 

G(J)/A(K)  >  0 

^(J)“P(K)*A(K)<0  n 

Y  G(J)-P(K)*A(K)<0/N 

Dl:- 

(C(K)-F(K))*A(K)/ 

(G(J)-P(K)*A(K)) 

Dl:- 

(C(R)-F(K) )*A(K) / 
(G(J)  -P(K)*A(K) 

D2:-  -F(K)*A(K)/G(J) 

5  TART 

- > COMP ARE 

Subroutine  MFL0G2 
Figure  5.7 


5.  (MIRROR)  Sec  j  Co 
Che  origin  of  node 
-k.  CalculaCe  D1 
and  D2  as  appro - 
prlaCe.  Go  Co  6. 

6.  (COMPARE)  If  Che 
new  value  of  D1  or 
D2  Is  less  Chan  MF, 
replace  MF  wich  Che 
new  value  and  keep 
crack  of  Che  leav¬ 
ing  arc.  Go  co  2. 


MIRROR 


J:-0(-R) 


G(J)  <  0 


^(J)  +  P(-K)  < 

N  G(J)  +  P(-K) 

Dl:» 

Dl:- 

(F(-K)-C(-K))/ 

(F(-K)-C(-K))/ 

(G(J)+P(-K)) 

(G(J)+P(-K)) 

D2:-F(-K)/G(J) 

ii 


->COMPARE 


Subrouclne  MFL0G2 
Figure  5.7(conC) 
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If  we  consider  arc  11,  we  are  in  SUBCASE  A  of  case  1  where  k  >  0 
and  Yj  <  0.  We  note  that  ^ j  -  p^a^  <_  0  and  therefore 

rn  "  _fiiaii/  y7 

-  1.2. 

We  may  determine  the  remaining*  reference  numbers  for  the  basic  arcs 
by  the  methods  discussed  in  this  section  and  these  values  are: 

r5  -  3.941 
r_6  -  0.9529 
r7  -  4.608 
r_9  -  0.8257 
r12  -  0.1482 
r_19-  0.5482 

The  minimum  occurs  for  arc  12.  Therefore  the  maximum  reference 
number  for  the  basic  arcs  is 

r  -  0.1482 


5.7  Determining  the  Entering  Arc 

In  the  previous  sections  we  developed  the  mechanisms  to  find 
an  arc  to  leave  the  basis.  Since  this  is  a  dual  algorithm,  we  now 
need  to  develop  the  techniques  to  be  able  to  select  an  admissible 
entering  arc. 

5.7.1  Characteristics  of  the  Entering  Arc 

An  entering  arc,  kg,  must  have  two  characteristics  to  enable 
it  to  enter  the  network.  The  first  characteristic  is  that  the 
entering  arc  must  be  admissible.  Since  the  algorithm  is  applied  to 
the  marginal  network  representation  that  admits  both  forward  (k  >  0) 
and  mirror  (k  <  0)  arcs,  we  define  the  following  admissibility  rule: 

k  >  0  and  c^  -  f^  >  0  then  arc  k  is  admissible 
k  <  0  and  f_fc  >  0  then  arc  -k  is 

admissible. 

We  can  state  this  in  a  more  intuitive  manner  as 

a  forward  arc  is  admissible  if  the  flow  in  the 


original  arc  is  less  than  its  capacity  and  a 


mirror  arc  Is  admissible  If  the  flow  in  the 

original  arc  is  greater  than  zero. 

The  second  characteristic  of  the  entering  arc  concerns  its 
location  in  the  network  relative  to  the  leaving  arc.  When  the 
leaving  arc  is  deleted  from  the  basis,  it  divides  the  nodes  of  the 
network  into  sets  and  N2>  The  set  of  nodes  which  remains 
connected  to  the  terminal  node  of  the  leaving  arc  k^  (j^)  is 
designated  the  N2  set.  We  will  denote  the  basis  subnetwork 
containing  as  D2«  We  note  that  always  forms  a  tree  rooted  at 
the  terminal  node  of  the  leaving  arc  and  no  longer  includes  the 
slack  node  or  a  cycle.  The  set  is  the  complementary  set  of 
nodes,  that  is  ■  N  -  N2.  We  denote  the  basis  subnetwork  that 
contains  as  .  Thus  consists  of  one  or  more  components  that 
always  includes  the  slack  node  and  any  flow-generating  cycles  that 
remain  after  we  delete  the  leaving  arc  from  the  basis. 

We  may  mathematically  express  the  components  that  are  formed 
when  *s  removed  from  the  optimal  basis  tree,  D  *  [N.M^,]  as 

Dl  -  [Nj.Mj] 

d2  -  [n2,m2] 

where 


iLe  N1  hz  N2 


N  U  N2  -  N 
Mi  U  »  Hj,  ~ 


When  arc  k^,  which  may  be  a  forward  arc  (k^  >  0)  or  a  mirror 
arc  (k^  <  0)  depending  on  Che  orientation  required  for  the  spanning 
tree,  is  deleted  from  the  basis,  there  no  longer  exists  a  way  to 
continue  augmenting  flow  on  the  candidate  arcs.  Thus  the  problem  is 
to  find  an  admissible  arc  to  add  to  the  basis  network  that  will 
reestablish  a  flow  augmentation  path  to  node  jL  and  therefore  allow 
continued  flow  augmentation  on  the  candidate  arcs.  If  we  let  k£  M^, 
where  is  the  set  of  admissible  arcs,  k  may  enter  the  basis  in  one 
of  three  ways.  First  k  may  join  Ni  to  ^  which  we  may  write  as  ke 
(NpN2).  The  notation  (NpNj,)  defines  the  set  of  arcs  each 
originating  at  a  node  in  the  set  Ni  and  terminating  at  node  in  N2« 
Secondly,  k  may  join  N2  to  Np  which  we  write  as  ke(N2,Ni).  Lastly 
k  could  join  N2  to  N2,  which  we  write  as  ke(N2,N2>.  In  the  first 
two  cases  the  subnetwork  D2  is  joined  to  a  tree  rooted  at  the  slack 
node  or  at  a  cycle  in  Dp  In  the  later  case  a  new  component  of  the 
basis  network  is  formed  from  D2  and  kg.  This  component  is  rooted  at 
a  new  cycle. 


We  now  can  describe  the  second  characteristic  that 
determines  the  selection  of  the  entering  arc.  It  oust  be  one  of  the 
following  four  cases,  which  is  illustrated  by  Figure  S.8. 

CASE  1:  1^  is  a  forward  arc  and  leaves  the  basis  at  its  upper  bound 

As  k^  goes  to  its  upper  bound,  we  may  think  of  additional 
flow  being  put  into  N2.  Thus  as  k^  leaves  the  basis,  k  must  form  a 
new  flow  augmenting  path  to  continue  to  allow  flow  to  be  put  into 
the  set  N2»  Thus  we  must  have  ke  (Nj,N2)  or  ke(N2,N2).  When  ke 
(N2,n2),  k  forms  a  flow  generating  cycle  (  0  >  1,  where  g  Is  the 
cycle  gain) . 

CASE  2:  k^  is  a  mirror  arc  and  leaves  the  basis  at  its  upper  bound 
As  k^  leaves  the  basis  by  going  to  its  upper  bound,  we  may 
think  of  flow  being  put  into  the  set  .  Thus  we  need  a  ke  (Hj.Hj) 
or  k£(N2,N2),  where  ke  (N2,N2)  forms  a  flow  absorbing  cycle  (6  < 

1,  where  8  is  the  cycle  gain). 

CASE  3:  k^  is  a  forward  arc  and  leaves  the  basis  at  its  lower  bound 
As  k^  leaves  the  basis  at  its  lower  bound,  we  are 
conceptually  putting  flow  into  the  set  N^,  Thus  we  must  have  an  arc 
ke  (N2,Nj)  or  k£  (N2,N2),  where  in  the  later  case  k  forms  a  flow 


absorbing  cycle.  (8  <  1,  where  8  is  the  cycle  gain). 

CASE  4:  k^  is  a  mirror  arc  and  leaves  the  basis  at  its  upper  bound 
We  can  again  think  of  flow  being  transferred  into  the  set 
leaves  the  basis.  Here  we  must  have  k£  (N^.N^)  or  k£  C^,^), 
and  in  the  later  case  k  forms  a  flow  generating  cycle  (  B>  1,  where  6 
is  the  cycle  gain).  This  will  allow  a  new  flow  augmentation  path  to 
put  new  flow  into  the  set. 

Cases  one  and  four  are  identical  as  to  the  form  of  ke  M. , 

A 

and  cases  two  and  three  are  the  same.  Table  S.l  summarizes  these 
results. 

5.7.2  On  Determining  the  New  Value  of  Y 

When  we  are  trying  to  determine  the  entering  arc,  we  must 
again  use  the  concept  of  node  gain.  However,  in  this  situation  we 
need  to  define  node  gain  differently  than  for  the  case  of  the 
leaving  arc.  In  this  context,  Yy  is  only  defined  for  the  nodes  in 
the  set  ^2’  In  this  case  Yu  is  defined  to  be  the  product  of  the 
arc  gains  starting  at  the  terminal  node  of  k^  (JL)  to  node  u,  where 
we  define  Y^  «  l  for  u  »  j^.  Thus  we  see  that  if  one  unit  of  flow 
arrives  at  and  travels  through  D2  to  node  u  the  flow  arriving  at 


node  u  is  Y  .  The  node  gains  can  be  computed  from  the  structure  of 

the  basis  tree  and  the  value  of  Y  .  For  any  arc  k(i,j)  that  is 

JL. 

included  in  the  tree  we  have 

Y.  -  Y^  for  ke  M2  (5.15) 

i  ■  o(k) 

j  -  t(k) . 

Node  JL  is  the  root  of  a  tree  in  Nj.  Assume  we  have  listed  the 
arcs  of  M2  in  predecessor  order  (  that  is  an  arc  appears  in  the  list 
after  all  its  predecessors).  Since  we  are  given  the  value  of  Y 

j  i* 

*  1,  progressing  through  the  ordered  list  Mj,  we  can  use  equation 
(5.15)  to  compute  the  node  gains  for  each  member  of  N2»  ROOTG 
(46 .page  247]  obtains  the  ordered  list  M2  and  NEWGAM,  as  part  of 
FIND2G,  calculates  the  Y  values. 

We  note  that  if  k£  (N2,N2),  arc  k  will  always  form  a  cycle 
with  gain  of 

6  -  ^  Yj/  Yt.  (5.16) 

This  equation  will  be  useful  in  the  algorithms  described  because  it 
allows  the  calculation  of  the  cycle  gain  without  an  expllct 
identification  of  the  arcs  on  a  cycle. 


As  an  example  consider  Figure  5.1(b)  in  which  arc  2(1,3) 


leaves  the  basis.  Here 


are 


3  and  the  value  of  the  node  gains  in 


i 

i 


I 

I 


Y3a7  "  a7 
Y2a5  *  a7a5 
Y5all“  a7a5al 1 
Y  7a5  «  a7a5al ia15 
Y5a-17  *  a7a5/a17 
Y3a6  *  a6 
V-9  *  a6/a9 
Vl4  *  a6al4* 


5.1.3  Selecting  the  Entering  Arc 


In  order  to  maintain  the  optimality  of  the  basis,  we  must 
select  an  arc  k  to  enter  the  basis  which  satisfies  the  conditions  of 


5.7.1  and  has  the  smallest  marginal  change  in  cost  per  unit  flow 
change  in  the  arc  k.  In  order  to  do  this,  we  must  develop  the 
marginal  cost  equations.  We  will  first  develop  the  equations  which 
deal  with  k  e  (NpNj)  or  k  £  (Nj.N^).  After  developing  these 
formulas,  we  will  derive  the  formulas  to  calculate  the  marginal  cost 


t 


I 
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for  k  e  (Nj,^) . 

CASE  1:  k^  is  a  forward  arc  and  is  going  to  an  upper  bound  or  k^  is 
a  mirror  arc  and  is  going  to  a  lower  bound  and  ke  MA  is  a  forward 
arc. 


If  we  let  d^  be  the  change  in  cost  per  unit  flow  change  in 
the  arc  k(i,j)  where  ketNj,^),  as  illustrated  in  Figure  5.9,  we 
can  calculate  d^  as 

dk  *  [(  ^i  +  \)/ak  _7r  j  1  Y  j*  (5.17) 

This  formula  is  used  in  the  subroutines  ABSORB  and  GENER2. 

CASE  2:  k^  is  a  forward  arc  and  is  going  to  an  upper  bound  or  k^  is 
a  mirror  arc  and  is  going  to  a  lower  bound  and  ke  MA  is  a  mirror 
arc. 

For  arc  ke  MA  and  k  <  0  where  k  e  (Nj.Nj),  as  illustrated  in 
Figure  5.10,  we  obtain 


[<*j  +  h-k)/a-k  "  *i,Y  i 


Recalling  h_k  *  -h^/ and  a 


l/-a^  in  an  expanded  network. 


d-k  -  (’A  -  "k  -’i)Y  i- 


(5.18) 


Equation  (5.18)  is  used  in  the  subroutines  ABSORB  and  GENER2. 


CASE  3:  k^  is  a  forward  arc  and  is  going  to  a  lower  bound  or  k^  is 
a  mirror  arc  and  is  going  to  an  upper  bound  and  ke  is  a  forward 
arc. 

The  marginal  cost  for  ke(N2»Nj)  per  unit  flow  change  in  arc 

k  is 


dk'l(,i*hk'*k,,jI  i- Vk1 
•  <  "A  '  "k  "*  1>  Y1 

which  is  Identical  to  equation  (5.18).  This  is  illustrated  in  Figure 
5.11. 

CASE  4:  k^  is  a  forward  arc  and  is  going  to  a  lower  bound  or  k^  is 
a  mirror  arc  and  is  going  to  an  upper  bound  and  k£  is  a  mirror 
arc. 

The  formula  for  the  marginal  cost  for  ke  (N^Nj),  k  <  0,  as 


ilustrated  in  Figure  5.12,  is 


which  is  Che  same  as  equation  (5.17). 

CASE  5:  Arc  k  e  Ma  and  ke(N2,N2). 

In  this  case  an  admissible  arc  k(i,j)  where  k£  (N2,N2)  forms 
a  cycle  as  shown  in  Figure  5.13.  tfe  must  insure  the  new  cycle  has 
a  gain  different  than  one.  Next  we  must  evaluate  the  marginal  cost 
of  the  path  formed  with  the  addition  of  arc  k.  The  path  forms  four 
parts. 

P,  :  the  path  from  some  junction  node  1  to  i 
P2  :  the  added  arc  k 

Pj  :  the  path  from  j  to  the  junction  node  1 
P^  :  the  path  from  j^  to  1. 

If  one  unit  of  flow  starts  at  node  j^»  this  will  require  a  flow  of 
^  at  node  1.  The  flow  into  node  1  from  the  circuit  P^  must  be  [6  / 
(6  -  1) ]  Y  ^  and  the  flow  entering  P^  from  P^  and  P^  at  node  1  must 
be  /  (  6-  1),  by  conservation  of  flow  arguments.  The  flow  out  of 
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-  1)1  *  [  Yj  /  ].  The  cost  of  the  four  parts  is: 


P3  -  [  8/  (6-  01  [  7T1  Y1  -TT  i  Yj  ] 

P2  -  t  0/  (  B-  O]  [  Yj>  /  ajc  1 
Px  *[8/(8-  Ol  [  Tr£  ~  (akir  x  Yj_)/  8  Yj  ]  Y^/a^ 

Adding  these  four  parts  and  simplifying  allows  us  to  write 

dk  -  [  8  /  (  8-  Ol  [(  *1  +  hk)/ak  _7r  -j  J  *Y  j  (5'19) 

where  8-  (^Y^/Y^  for  ke(N2,N2).  This  is  the  formula  used  in 
ABSORB  and  GENGR2. 

For  admissible  mirror  arcs,  k(j,i)  k  <  0,  the  value  of  d_fc 

becomes 

d_k  -  (6  /(  8-  01  [  ff  jak  -hk-^1lYi  (5.20) 

where  8-  Y^  /  (^  Y^  )  f0r  -ke(N2,N2).  This  is  the  formula  used 
in  ABSORB  and  GENER2. 


Recall  that  in  section  5.7.1  we  outlined  four  cases  which 
determined  the  characteristics  of  the  entering  arc.  In  cases  1  and 
4  we  needed  a  flow  generating  cycle.  Using  equation  (5.16)  we  may 
determine  if  arc  k  has  g  >  1,  g<  1,  or  8*1.  If  8*1,  arc  k  is  not 
admissible.  If  arc  k  has  g  <1.  and  the  mirror  arc  k,  k  <  0,  is 
admissible,  we  will  calculate  the  marginal  cost  for  arc  k.  If  arc  k 
has  g  >  1  we  calculate  the  marginal  cost  for  arc  k,  k  >  0.  The 
similar  procedure  is  used  for  cases  2  and  3  where  the  entering  arc 
needs  to  form  a  flow  absorbing  cycle. 

Since  d^  represents  the  marginal  increase  in  cost  per  unit 
flow  change  in  an  admissible  arc,  we  need  a  dual  type  iteration  to 
find  the  arc  with  the  smallest  value  of  d^.  This  will  be  the 
entering  arc  and  receive  the  designation  kg.  Recalling  if  is 
the  set  of  admissible  arcs  from  the  expanded  network  such  that 

if  k(i,j)  e  either 

ke(N1,N2),  ke(N2,N1),  or  ke  (N2,N2). 

Then  the  entering  arc  is  determined  by 


PURPOSE:  To  find  (  PIND2G  ) 

and  arc  to  enter 

the  basis.  IHITIALl 


1  .(INITIAL 1)  Set 
node  set  S(I)-O.If 
KL  is  a  mirror  arc 
find  tree  rooted  at 
O(-KL).  If  not  find 
tree  rooted  at 
T(KL) .  If  tree 
rooted  at  IROOT 
forms  a  cycle,  de¬ 
lete  the  second  re¬ 
ference  to  node  j 
from  LN  and  delete 
KL  from  LA. 

• 

2 .  ( INITIAL2 )  For 
each  node  in  N2, 
set  S(I)  -  1.  In¬ 
itialize  DL  and 
KE. 

3. (NEWGAM)  Set  T  for 
root  node  equal  1. 
Compute  and  assign 

7*8  for  each  node 
in  N2. 


FOR  I:»l  TO  H 


4.  (SEARCH)  Start 
search  through  list 
of  arcs.  If  k  is  a 
basic  arc  or  an 
artificial  arc  go 
to  8.  Otherwise  go 
to  5. 

5.  (DECIDE)  Set  I 
equal  to  the  origin 
node  of  arc  k.  If 
k  is  a  forward  arc 
go  to  6,  otherwise 
go  to  7. 


INITIAL2 


S(IRQOT) :-l 


ic  -  o 


FOR  L:-l  TO  IC 


Ijj:-  LN(L+l),  S(JJ):-1 


DL:-  999999,  KE:-  0 


- >  NEWGAM2 


Subroutine  FIND2G 
Figure  5.14 
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6.  (FORWARD)  If  KL 
is  stored  as  a  for¬ 
ward  arc  and  KL  Is 
going  to  a  lower 
bound  (g(j)  <  0), 
call  GENER2 .  Other¬ 
wise  call  ABSORB. 

Go  to  8. 

7.  (MIRROR)  If  KL  is 
stored  as  a  mirror 
arc  and  going  to  an 
upper  bound  (g(i)  < 
0),  call  GENER2. 
Otherwise  call  AB¬ 
SORB.  Go  to  8. 

8.  (COUNT)  When  all 
the  arcs  have  been 
checked,  call 
PIVOT1G  to  pivot 

in  KE,  otherwise  go 
to  A. 


NEWGAM2 


G(IR00T):-1 

Y\ 

IC  -  0 

/* 

// 

// 

/  / 

FOR  L:-l  TO  IC 

K:«LA(L),  JJ:-LN(L+1) 

/  / 

// 

/  / 

Y\^  K 

>  0  ^/N 

// 

// 

// 

II :*0(K) 

G(JJ):«G(II)*A(K) 

II:-T(-K) 

G(JJ):-G(II)/A(-K) 

- >SEARCH 

SEARCH 


Subroutine  FIND2G 
Figure  5.14(cont) 
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MIRROR 


G(KI)>0 


C  ABSORB  ) 


(  GEMER2  ) 


— >COUNT 


COUNT 


Subroutine  FIND2G 
Figure  5»l4(cont) 


1.  (SEARCH)  If  arc 
k  originates  In  N1 
( S ( I ) “0 )  and  term¬ 
inates  In  N2  (S(I)- 
1),  check  forward 
arc  for  admissibil¬ 
ity.  If  admissible 
go  to  2.  If  arc  k 
originates  in  N2 
and  terminates  in 
Nl,  check  mirror 
arc  for  admissibil¬ 
ity.  If  admissible 
go  to  3,  otherwise 
go  to  4. 

2.  (FORWARD)  Evalu¬ 
ate  dk  for  forward 
arc  and  go  to  7. 

3.  (MIRROR)  Evaluate 
dk  for  mirror  arc 
and  go  to  7. 

4. (CYCLE1)  Compute 
gain  of  cycle  by 
inserting  forward 
arc. If  0  »  1,  re¬ 
turn.  If  gain  >  1, 
test  forward  arc 
for  admissibility. 
If  admissible ,go  to 

5.  otherwise  re¬ 
turn.  Test  mirror 
arc  for  admissi¬ 
bility.  If  admis- 


MIRROR 


>:-(PI(J)*A(K)-H(K)-PI(I))*G(I),  KK:— K 
- COMPARE 


CYCLE 1 


BBT:-A(K)*G(I)/G(J) 


BET-1 


BET  >  1 


F(K)  <  C(K) 


Y\  P(K)  >  0 
BET:-i/BET| 


RETURN  TWO CYC  RETURN  MIRCYC  RETURN 


Subroutine  ABSORB 
Figure  5.15 


sible,  go  co  6, 
otherwise  recurn. 

5.  (FWDCYC)  Coapuce 
dk  for  forward  arc 
chac  foraa  a  cycle. 
Go  to  7. 

6. (MIRCYC)  Coapuce 
dk  for  airror  arc 
chac  foraa  a  cycle. 
Go  Co  7. 

7.  (COMPARE)  Test  dk 
againsc  saallest 
found  co  chi a 
point.  If  dk  ia 
saaller,  call  chia 
che  entering  arc. 


FWDCYC 


D:«(BET/ (BET-i))*(((PI(I)+H(K)/A(K))-PI(J))* 
G(J) 


IKK:-  K 


- COMPARE 


MIRCYC 


D:-(BET/(BET-1))*(PI(J)*A(K)-H(K)-PI(J))*G(I) 


KK:— K 


! 


COMPARE 


3 


PURPOSE:  To  deter¬ 
mine  if  are  k  is 
admissible  to  en¬ 
ter  the  basis  and 
calculate  the 
minimum  dk  value. 

1.  (SEARCH)  Check 
if  arc  k  origi¬ 
nates  in  N1  (S(I) 
•0)  and  terminates 
in  N2  (S(I)-i).  If 
so,  check  mirror 
arc  for  admissibil¬ 
ity.  If  admissible 
go  to  2,  otherwise 
return.  If  arc  k 
originates  in  N2 
and  terminates  in 
Nl,  check  forward 
arc  for  admissi¬ 
bility.  If  ad¬ 
missible  go  to  3, 
otherwise  return. 

If  arc  k  goes  from 
N2  to  N2,  go  to  4. 

2.  (FORWARD)  Evalu¬ 
ate  dk  for  forward 
arc.  Go  to  7. 

3.  (MIRROR)  Evalu¬ 
ate  dk  for  mirror 
arc.  Go  to  7. 

4. (CYCLE2)  Com¬ 
pute  gain  caused 
by  adding  arc  k.  If 
cycle  gain  >  1  and 
mirror  are  is  ad¬ 
missible,  go  to  6. 
If  cycle  gain  is 

<  1  check  forward 
arc  for  admissib¬ 
ility.  If  admissi¬ 
ble  go  to  5. 
Otherwise  return. 


SEARCH 


S(I)-0 


\ 


\  S(J)-0 


TURN  RETURN 


FORWARD 


>:-(((PI(I)+H(K))*A(K)-PI(J))*G(J) 


IKK:-  K 


- >C0MPARE 


MIRROR 


.(PI(J)*A(K)-H(K)-PI(I))*G(I) 


C:-  -K 


- COMPARE 


CYCLE 2 


RETURN 


BET : -A(K) *G( I) /G( J ) 


BET  -  1 


£ 


BET  >  1 


YV  F(K)  >  0 
BET:-1/BET 


,F(K)  <  C(K)y 


RETURN  MIRCYC  RETURN  FVDCYC  RETURN 

Subroutine  GENER2 
Figure  5.16 


I 


I 


I 

s 
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dk  -  Min  {dk}.  (5.21) 

k£  Ma 


This  is  the  formula  use  in  the  subroutine  COMPARE.  It  is  clear  that 

d,  is  non-negative  because  the  process  begins  with  a  dual  feasible 
KE 

solution. 

The  concepts  developed  above  are  implemented  in  the 
subroutines  FIND2G,  ABSORB,  and  GENER2. 


5 ,8  A  Basis  Change  Method 


|  At  each  iteration  subroutine  MFL0G2.  will  determine  the  arc  ~~~* 

that  limits  the  flow  change  on  the  flow  augmenting  trail(s),  and 

that  arc  will  leave  the  basis.  Subroutine  FIND2G  will  determine  the 

|  entering  arc.  Arc  is  such  that  k^  e  (N^N^),  kgE  (N2,N2),  or  kgE 

(N2,Nj).  In  the  first  two  cases  the  subroutine  TRECHG  [pg  116,(46)] 

can  be  used  directly  to  accomplish  the  basis  change  operation  for 

-  the  PSENG  algorithm.  However,  if  k^e  (N2,N^) ,  kg  must  enter  as  -k£ 

as  TRECHG  assumes  k^'s  origin  is  in  the  set  N^. 

In  order  to  update  the  dual  variables,  tt  ,  the  set  of  nodes 

and  arcs  in  the  tree  rooted  at  j,  is  found  with  the  subroutine 

■  *E  - - 

ROOTG  [pg.247,  (46)]. 

I 

| 

. _ '  -  _ _ •  . _ - _ _ .  -  - _ . _ _ -  ■ _ ' _ _  a _ _ _ -  _ ..a.  a  - _ aa _ a  ■  ,  -  - _  -  -  ■  -  - _ -  _ _ _ -  ■  - _ 1 _ a _ a _ - _ -  _ a _ a. _ a  - ■  -a  -  a- - a _ a  a  M 


PIV0T1G 


PURPOSE:  To  change 
Che  basis  tree  by 
adding  KE  and  de¬ 
leting  KL.  Com¬ 
pute  new  w  values 
for  new  basis. 


Subroutine  PIV0T1G 
Figure  5.17 


If  kg  has  formed  a  cycle,  subroutine  CYCLE  [ pg .273,  (46)]  is 
used  to  calculate  the  gain  of  the  cycle  formed  by  k^.  The  n  values 
are  then  updated  with  the  subroutine  DUAL  [pg.271,  (46)].  These 
concepts  are  implemented  in  the  subroutine  PIV0T1G. 


5.9  The  Complete  Algorithm 

At  this  point  we  can  state  the  complete  algorithm  used  in  PSENG. 

1.  Solve  the  original  problem  to  optimality  by 
any  primal  or  dual  solution  technique. 

2.  Determine  the  list  of  candidate  arcs  and 
the  list  of  parametric  parameters. 

3.  Find  the  flow  augmenting  trail  or  trails 
which  are  generated  by  the  candidate  arcs. 

Determine  the  arc  to  leave  the  basis.  If 
none  exists,  go  to  step  4.  If  the  leaving 
arc  is  a  non-basic  candidate  arc,  go  to 
step  5.  Otherwise  go  to  step  6. 

4.  At  least  one  basic  candidate  arc  can  have 


its  capacity  increased  without  bound 
indicating  no  further  parametric  changes 


will  affect  the  basis.  Stop. 

5.  Augment  the  flow  as  much  as  possible  on  the 
trail  or  trails.  Adjust  the  capacity  and 
flows  on  the  candidate  arcs.  Stop,  any 
futher  iterations  would  cause 
infeasibility. 

6.  Augment  the  flow  as  much  as  possible 
through  the  trail  or  trails.  Adjust  the 
capacity  and  flows  on  the  candidate  arcs. 

7.  Find  an  arc  to  enter  the  basis. 

8.  Change  the  basis  by  deleting  the  leaving 
arc  and  inserting  the  entering  arc.  Modify 
the  dual  variables  to  satisfy  complementary 
slackness.  Return  to  step  3. 

It  is  important  to  remember  that  for  any  gains  algorithm 
that  the  arithmetic  is  performed  with  real  numbers  rather  than  with 
Integers.  Thus  the  problem  of  computer  round-off  error  is  quite 
likely  to  occur.  The  tests  for  equality  must  take  this  into 
account.  For  example  in  the  subroutine  ABSORB,  the  test  f^  >  0  is 
replaced  by  f^  >  e  , where  e  is  a  small  quantity.  While  the  flow 
charts  presented  in  this  section  and  in  previous  sections  do  not 
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4. (CHECK)  Initial¬ 
ize  T's  to  zero 
and  node  check  val¬ 
ues  to  zero.  For 
each  candidate  arc 
check  if  the  arc  is 
nonbaslc  and  its 
parametric  para¬ 
meter  is  negative. 
If  so  calculate  the 
reference  number. 
Store  smallest  re¬ 
ference  number  as 
HE.  If  parametric 
parameter  is  neg¬ 
ative  set  IE  to  the 
terminal  node  of  XK 
otherwise  let  IE  be 
the  origin  node. 
Call  PATH  which 
will  create  LA,LN, 
and  RC  lists  for 
all  flow  augment¬ 
ing  trails  created 
by  the  nonbaslc 
candidate  arcs.  If 
no  cycles  are  found 
(z-0)  go  to  5. 
Otherwise  call 
CYCLEG  to  update 
T’s  for  affects 
of  the  cycle. 


CHECK 


Main  Program  PSENG 
Figure  5.18(cont) 


MFL0G2  to  deter¬ 
mine  the  maximum 
reference  number 
for  the  basic  arcs. 
If  no  basic  arc  can 
be  found  to  leave 
the  basis,  call 
COMPLETE  and  re¬ 
turn  to  MAIN  MENU. 
Otherwise  call  UP¬ 
DATE  to  change  arc 
flows  and  candidate 
arc  capacities.  If 
in  update  a  stop¬ 
ping  condition  is 
reached  (st»l), 
call  RESET  and  re¬ 
turn  to  MAIN  MENU. 
Otherwise  go  to  4. 


ITER: “ITER+1 ,  MF : -ME 

(  MFL0G2(IC,KL,MF)  ) 

KL  -  0  yAl 

(COMPLETE  (MF.  I.  XM)) 

(  UPDATE (MF,KL)  ) 

(  RESET  ) 

ST  -  1  yAi. 

(  RESET  ) 

- >  MAIN  MENU 

- >MAIN  MENU 

Main  Program  PSENG 
Figure  5.18(cont) 
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PURPOSE:  To  update 
flows  and  capacit¬ 
ies  on  candidate 
arcs  and  to  update 
flows  on  appropri¬ 
ate  basic  arcs* 

Also  to  find  KE  to 
enter  basis,  and 
to  print  out  in¬ 
termediate  results. 

If  MF  ■  0,  this  is 
a  degenerate  iter¬ 
ation.  Otherwise 
update  flows  for 
the  basic  arcs 
found  in  LA.  Up¬ 
date  the  f lows  and 
capacities  for  the 
candidate  arcs* 

Then  print  itera¬ 
tion  results.  Call 
FIND2G  to  determine 
KE,  the  entering 
arc.  If  a  stopping 
condition  has  been 
found  call  PRINT2 
to  print  final  re- 
[  suits ,  otherwise 

return. 


C  UPDATE (MF,KL)  ) 


(  FLOWGCLA.LN, IC.MF)  ) 

/  / 

FOR  I:-l  TO  XM  ' 

// 

// 

// 

// 

// 

// 

// 

// 

XK:-XK(I),ZI:-0(XK),EJ:-T(XK), 

C(XK) : -C(XK)-rtff*PP(I) 

Vy  PB(ZJ)«XK  OR  PB(ZI)— XK  yZ 

F(XK)-0  AND  P(XK)  >  0 ^/N 

F( XK) : “F(XK)+MF*PP ( I ) 

(  PRINT1  ( XK( I ) , PP ( I ) , CL( I ) , C ( I ) , F ( I ) , H( I ) ) 


(  FIND2G  ) 


[ PRINT2(I.AH(I),XL(I),MF(I))  ) 


RETURN 


Subroutine  UPDATE 
Figure  5.19 


E 


B 

I 

G 


& 


PURPOSE:  Final  up¬ 
date  of  flows  on 
basic  arcs  and  to 
update  capacities 
on  candidate  arcs. 
Print  final  iter¬ 
ation  results  and 
stannary  table. 


Note:  If  MF  • 
99999,  a  basic  arc 
can  have  its  cap¬ 
acity  increased 
without  bound. 


C  COMPLETE (MF. I. XM)  ) 


(  FLOWG(LA.LN, IC .KL.ILC.MF)  ) 


// 

// 

// 

// 

// 

// 

FOR  I:-l  TO  XM 

XK: -XK( I ) , 21 : -O(XR) , ZJ : -T(XK) 
C(XK):-C(XK)  +  MF*PP(I) 

K  PBfZJ)-XK  OR  PBfZI)—  XK 

F(XK):-F(XK)  +  MF*PP(I) 

o 

a 

KL:«NK 


C  PRINTlCXK,PP(I),CLCI),C(I)tF(I)tH(iyD 
C PRIMT2(I,Att(l) .XL(I).MF(I) )  ) 

RETURN 


Subroutine  COMPLETE 
Figure  5.20 


reflect  such  adjustments,  these  adjustments  do  appear  In  the 
computer  implementation  of  PSENG.  The  concepts  presented  above  are 
implemented  in  the  algorithm  PSENG. 

Figure  5.21  presents  an  example  application  of  the  algorithm 
PSENG.  Figure  5.21(a)  lists  all  the  initial  parameters  of  the 
problem  and  Figure  5.21(b)  shows  the  optimal  basic  feasible  solution 
derived  by  a  primal  solution  method.  In  each  of  the  Figures  5.21(c) 
through  5.21(h),  Che  entering  arc  is  indicated  by  a  dashed  line  and 
the  leaving  arc  is  marked  by  an  X.  The  basic  arcs  are  indicated  by 
heavy  lines  and  the  non-basic  arcs  at  capacity  by  light  lines.  The 
dashed-uct  lines  indicate  non-basic  candidate  arcs  with  flows  at 
capacity.  The  flow  and  capacity  at  the  completion  of  each  basis 
change  and  the  resulting  total  cost  of  those  flows  are  also  given  at 
each  step. 

The  candidate  arc  list  and  their  respective  parametric 
parameters  are  listed  below 

X  -  [5,6,7,11,14,15,20] 

P  -  (-.3, 1.2, -2.1, -.6667, .1,-3, 1.4]. 

In  Figure  5.21(c)  candidate  arc  20  creates  two  flow  augmenting 
paths.  The  values  of  theY's  are  computed  in  PATH,  and  MFL0WG2 
calculates  a  reference  number  of  1.10.  Arc  2  is  forced  to  leave  the 


2.20.1) 


m 


(ck,hk>ak) 

(bj,b$j,hS{  ) 


<fk-V 

cost  =132.400 


8,10) 


(Vk) 


cost  =  124.060 


^k,ck  ) 


cost  -  129.009 


basis  at  its  upper  bound.  The  subroutine  UPDATE  calls  FIND2G  which 
determines  the  entering  arc,  k^  -  3,  and  updates  the  flows  and 
capacities  on  the  candidate  arcs  and  the  flows  on  the  basic  arcs  on 
the  flow  augmenting  trails. 

Figure  5.21(f)  shows  four  augmenting  trails:  two  generated 
by  non-basic  candidate  arc  14  and  two  generated  by  non-basic 
candidate  arc  15.  Figure  5.21(h)  shows  a  cycle  which  is  formed  when 
arc  12  enters  the  basis.  The  algorithm  terminates  when  the  flow  and 
capacity  of  candidate  arc  15  are  driven  to  zero.  This  is  shown  in 
Figure  5.21( j ) . 

Table  5.2  summarizes  the  results  of  the  algorithm.  Figure 
5.22  illustrates  that  a  parametric  analysis  on  the  capacities  of  a 
network  produces  a  piecewise-linear  and  convex  function.  It  is 
interesting  to  note  that  the  objective  function  decreases  until  the 
third  iteration,  when  the  objective  function  begins  to  Increase. 

5.10  Extensions  to  Node  External  Flows 

In  sections  5.2  through  5.9  we  examined  the  effects  of 
parametric  analysis  on  the  arc  capacity  vector.  In  the  context  of 
linear  programming,  we  have  only  analyzed  the  upper  bound  part  of 
the  right-hand  side  vector.  We  can  easily  extend  the  results  of 
sections  5.2  through  5.9  to  accommodate  the  effects  of  changing  the 
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ITERATION 

NUMBER 

OBJECTIVE 

FUNCTION 

VALUE 

LIMITING 

ARC 

0 

132.400 

1 

126.076 

02 

2 

124.339 

14 

3 

124.060 

15 

4 

128.009 

22 

5 

128.009 

17 

6 

133.606 

12 

7 

170.585 

18 

8 

180.188 

15 

Summary  Table  of  PSENG 
Table  5. 

Example  Problem 
2 

CUMULATIVE 

REFERENCE 

NUMBER 


1.096 

1.538 

2.530 

2.798 

2.798 

2.847 

3.792 

4.000 


supply  and  demand  values  expressed  by  Che  right-hand  side  vector  b. 
While  previously  we  focused  our  attention  on  changing  arc 
capacities,  we  will  now  show  how  to  extend  the  PATH  subroutine  to 
take  into  consideration  changes  in  the  external  flow. 

Because  PATH  starts  at  a  node  and  traces  the  flow  augmenting 
path  back  to  the  root  node  or  a  root  cycle,  we  merely  need  to 
designate  a  parametric  parameter  for  each  node  at  which  the  external 
flow  is  to  be  changed.  Then  PATH  will  determine  the  Y  's  for  each 
node  generated  by  these  flow  augmenting  trails.  Subroutine  MFL0WG2 
will  determine  the  limiting  arc.  We  then  need  only  determine  the 
entering  arc  in  FIND2G.  We  continue  this  process  until  any  further 
changes  will  cause  the  problem  to  become  infeasible  or  the  current 
basis  will  remain  optimal  for  any  further  changes. 


CHAPTER  6 


APPLICATIONS  OF  PARAMETRIC  SENSITIVITY  ANALYSIS 

6.1  Introduction 

In  this  chapter  we  develop  two  dual  Incremental  algorithms 
which  use  the  concepts  of  parametric  analysis  developed  In  Chapter 
5.  The  first  algorithm  Is  called  PARARC  which  uses  parametric 
analysis  on  arc  capacities.  The  second  algorithm  Is  called  PARANDE 
which  uses  parametric  analysis  on  node  external  flows.  The 
computational  efficiency  of  PARARC  and  PARANDE  are  compared  against 
a  dual-incremental  code  called  INCREMG  (46)  and  two  primal  codes, 
NETG  (38)  and  PGAINS  (46).  Lastly  we  discuss  applications  of 
parametric  analysis  to  other  classes  of  problems. 

6.2  Parametric  Arc  Incremental  Flow  Algorithm 

Given  a  network  with  arcs  having  capacity,  cost,  and  gain 
parameters  with  source  and  destination  nodes  defined  and  with  a 
required  quantity  of  flow  specified  at  the  sources  and  destinations, 
our  problem  is  to  determine  the  arc  flows  with  minimum  total  cost. 
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The  problem  can  be  stated  as 


Min  Z  h.  f,  (6.1) 

k=*l 

S.T.  -  Z  a^f^  +  £  ■  b.  ^or  a^’1  1  e  N“s  (6.2) 

keMn  keM0i 

0  £  fk  i  ck  ke  M  •  (6.3) 

The  parametric  arc  incremental  algorithm  is  a  dual-node  infeasible 
approach  applied  to  a  special  network.  The  network  has  a  single 
super-sink  (node  t)  and  a  single  super-source  (node  s). 

For  each  source  node  an  arc  is  created  from  the  super-source 
to  the  source  node  with  capacity  equal  to  b^  For  each  destination 
node  an  arc  is  creat  d  from  the  destination  node  to  the  super-sink 
with  capacity  equal  to  -b^,  We  let  the  expanded  arc  set  which  is 
created  be  denoted  as  mC  We  can  redefine  the  primal  problem  as 


Here  F  is  the  sum  of  the  destination  demands.  An  assignment  of 
flows  that  satisfy  (6.5)  and  (6.7)  will  be  called  F.  An  F  which 
satisfies  (6.4)  for  some  value  of  flow  at  the  sink  less  than  F  will 
be  called  an  intermediate  optimum.  The  flow  F  that  satisfies 
(6. 4-6. 7)  is  called  the  optimum. 

To  solve  this  problem  using  the  parametric  approach  of 
Chapter  5,  we  identify  the  arcs  from  the  destination  nodes  to  the 
super-sink  as  candidate  arcs  with  parametric  parameters  equal  to  the 
demand  at  the  destination.  Ue  then  increase  the  candidate  arc 
capacities  from  zero  until  they  reach  a  capacity  equal  to  the 
original  destination  demands.  This  creates  multiple  flow  augmenting 
paths  at  each  iteration  and  allows  flows  to  be  augmented  from  the 
super-source  to  the  super-sink. 

Thus  we  begin  with  an  intermediate  optimum  Fq  for  some  value 
of  output  flow  at  the  super-sink.  Next  we  obtain  a  shortest  path 
network  Dg*\  that  determines  for  each  node  the  minimum  cost  per  unit 
of  additional  flow  to  the  node  and  the  path  over  which  the  flow  may 
be  obtained.  The  output  flow  is  increased  in  D  ^  along  the  flow 
augmenting  trails  generated  by  the  candidate  arcs  until  one  or  more 
arcs  in  the  trails  become  capacitated.  This  flow  augments  Fq  to 
become  F^,  Then  F^  is  the  intermediate  optimum.  A  new  augmenting 
network  is  constructed,  Dg^,  and  the  output  flow  is  again  augmented 
to  obtain  F^.  The  process  continues  until  the  flow  at  the 
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super-sink  equals  Ft»  At  every  step  F^  is  an  intermediate  optimum; 
therefore  at  the  termination  the  flow  pattern  must  be  at  the 
optimum. 

The  general  steps  of  the  algorithm  are  as  follows: 


1.  Find  an  initial  least  cost  spanning  tree. 

2.  If  the  flow  at  the  super-sink  is  F£,  stop 
with  the  optimal  solution. 

3.  For  each  candidate  arc  find  the  flow 
augmenting  trails  which  Include  the  basic 
arcs.  Determine  the  maximum  allowable  flow 
Increment  which  will  occur  on  one  of  the 
trails.  Determine  the  leaving  arc. 

Augment  the  flows. 

4.  Find  an  arc  to  enter  the  basis.  If  none 
exists,  stop.  There  is  no  feasible 
solution.  Otherwise  go  to  5. 

5.  Change  the  basis  by  deleting  the  leaving 
arc  and  Inserting  the  entering  arc.  Modify 
the  dual  variables  to  satisfy  complementary 
slackness.  Return  to  step  2. 


We  now  will  develop  the  detailed  steps  of  the  algorithm. 


6.2.1  Obtaining  the  Initial  Basis  Tree 

In  order  to  apply  the  concepts  developed  in  Chapter  5,  we 
must  modify  how  the  network  is  to  be  represented.  Consider  the 
sample  network  problem  shown  in  Figure  6.1.  We  will  create  a 


super-sink  node,  which  will  have  a  fixed  external  flow  equal  to  the 
sum  of  the  destination  fixed  external  flows.  We  also  will  create 
"candidate  arcs”  from  the  destination(s)  to  the  super-sink  with  the 
following  parameters: 


(capacity) 


M  +  1 


(cost) 


■  -  {1^  l  o(k)  ]}  (parametric  parameter) 
a.  ■  1  (gain). 


We  will  also  create  arcs  from  the  slack  node  to  all  nodes  with  slack 
external  flows  with  the  arc  parameters 


ck  ’  bsi 

h.  -  0 
k 


pk  "  0  * 


Finally  we  create  an  artificial  arc  from  the  slack  node  to  the 
super-sink  with  parameters 
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hk“  M 

\m  1 
pk  -  0  . 

Since  PSENG  requires  conservation  of  flow  at  all  nodes,  an 
artificial  arc  is  required  to  obtain  conservation  of  flow  at  the 
super-sink.  When  the  flow  on  the  artificial  arc  goes  to  zero,  the 
flow  through  the  network  will  be  Fc.  All  of  these  transformations 
are  accomplished  in  the  subroutines  READG1  and  ORIGSGl. 

Once  the  transformed  network  is  obtained,  we  use  Dijkstra's 
algorithm  adapted  for  the  generalized  network,  as  documented  in 
[46, page  283],  to  obtain  an  initial  basic  feasible  solution.  The 
subroutine  TRE1NT  [46, page  121]  is  used  to  initialize  the  parameters 
of  the  basis  tree.  Since  a  generalized  network  flow  problem  may 
have  arc  gains  greater  than  one  or  have  negative  arc  costs, 
Dijkstra's  algorithm  may  not  provide  an  optimal  shortest  path  tree. 
Therefore,  we  now  use  the  algorithm  PSHRTG  [46, page  285]  to  find  the 
initial  least-cost  spanning  tree.  Once  this  tree  is  found  we  reset 
the  costs  on  the  "candidate  arcs”  to  be  zero.  This  is  illustrated 
in  Figure  6.3.  We  reset  these  costs  so  when  the  candidate  arc  flows 
Increase  from  zero  to  the  former  external  flows  of  the  original 
destination  nodes,  the  total  cost  for  the  solution  will  reflect  the 


actual  cost  for  the  problem.  At  this  point  we  have  an  initial 
optimal  (all  >_  0),  basic,  infeasible  (the  artificial  arc  has  flow 
>  0)  solution. 

6.2.2  Determining  the  Arc  to  Leave  the  Basis 

Similarly  as  in  PSENG,  we  now  use  the  PATH  subroutine  to 
determine  the  flow  augmenting  trails  for  each  candidate  arc.  This 
will  create  a  list  of  a; < 3  (LISA),  a  list  of  nodes  (LISN),  and  a 
list  of  gamma  values  (G).  Using  these  lists  MFLOG  [46, page  266] 
will  determine  the  maximum  flow  change  or  reference  number  on  these 
trails  and  the  limiting  arc.  We  will  denote  the  limiting  arc  as  k^. 

6.2.3  Determining  the  Entering  Arc 

We  use  the  subroutine  F1NDARC  to  determine  the  entering  arc, 
k^.  If  no  arc  can  be  found  to  enter  the  basis  at  this  iteration,  no 
feasible  solution  exists  to  the  proposed  problem. 

6.2.4  Changing  the  Basis 

Updating  of  the  basis  tree  representation  and  the  dual 


variables  is  accomplished  in  the  subroutine  PIV0T1G 


6.2.5  The  Complete  Algorithm 


We  can  now  give  a  complete  statement  of  the  parametric 
Incremental  algorithm. 

1.  Use  a  variant  of  Dijksta's  shortest  path 
algorithm  to  find  the  initial  least-cost 
spanning  tree.  Use  the  arcs  In  the 
spanning  tree  to  define  a  set  of  node 
labels  it  i  such  that  it  g  ■  0  and  tt  ^  «  (  it  ^ 

+  h^)  /  afe  for  all  arcs  k(i,j)e  Hj.. 

2.  Use  a  triple  labeling  scheme  to  label' and 
store  the  least-cost  spanning  tree. 

3.  If  the  flow  on  the  artificial  arc  from  the 
super-source  to  the  super-sink  node  is 
zero,  stop  with  the  optimal  solution. 

4.  For  each  candidate  arc  from  a  destination 
node  to  the  super-sink,  find  the  augmenting 


trails  which  include  the  candidate  arcs  and 
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Che  basic  arcs.  Determine  Che  maximum 
allowable  flow  increment  on  these  trails. 
Determine  the  arc  on  these  trails  which  is 
at  capacity  and  designate  it  the  leaving 
arc  k^.  Augment  the  flows. 

5.  Delete  arc  k^  from  the  spanning  tree, 

partitioning  the  network  into  two  subtrees 
defined  by  the  sets  and  Search  the 

set  of  admissible  arcs  originating  or 
terminating  in  which  connect  to  the  set 
N2  or  the  admissible  arcs  that  originate 
and  terminate  in  .  If  no  admissible  arc 
exists,  stop;  there  is  no  feasible 
solution.  Otherwise  select  the  arc  with 
the  minimum  d^  and  designate  it  as  the 
entering  arc  kg. 

6.  Perform  a  basis  change  which  determines  the 
new  least-cost  spanning  tree.  The  basis 
change  is  accomplished  by  deleting  k^  and 
adding  kg.  Update  the  dual  variables  and 
return  to  step  2. 
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The  algorithm  PARARC  implements  these  concepts.  The  flow  charts  for 
the  FORTRAN  coding  of  this  algorithm  can  be  found  in  Appendix  1. 

Figure  6. A  presents  an  example  application  of  the  algorithm 
PARARC.  Figure  6.4(a)  shows  the  initial  problem  data.  Figure 
6.4(b)  illustrates  the  intermediate  transformation  which  creates  the 
four  candidate  arcs  from  the  original  destinations  to  the  super-sink 
and  the  creation  of  the  artificial  arc  from  the  slack  node  to  the 
super-sink. 

In  each  of  the  Figures  6.4(c)  through  6.4(k)  the  entering 

* 

arc  is  indicated  by  a  dashed  line  and  the  leaving  arc  by  an  X.  The 
basic  arcs  are  indicated  by  heavy  lines  and  the  nonbasic  arcs  at 
capacity  by  light  lines.  The  dashed-dot  lines  indicate  the 
non-basic  candidate  arcs.  The  flow  and  capacity  at  the  completion 
of  each  basis  change  are  given  at  each  step. 

Figure  6.4(c)  shows  the  initial  least-cost  spanning  tree 
derived  from  a  variant  of  Dijstra's  algorithm.  Figures  6.4(d) 
through  6.4(j)  show  the  intermediate  steps  in  the  solution  process 
as  the  flow  on  the  artificial  arc  is  being  decreased.  It  is 
interesting  to  note  at  each  intermediate  iteration,  we  are 
increasing  the  flows  on  multiple  flow  augmenting  paths,  rather  than 
a  single  path  as  in  the  usual  flow  augmentation  approach  (46). 


Example  Application  of  PARARC  Algorithm 
Figure  6.4(cont) 
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Example  Application  of  PARARC  Algorithm 
Figure  6.4(cont) 


Example  Application  of  PARARC  Algorithm 
Figure  6.4(cont) 
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Example  Application  of  PARARC  Algorithm 
Figure  6.4(cont) 


Example  Application  of  PARARC  Algorithm 
Figure  6.4(cont) 


Figure  6.4(k)  shows  the  optimal  solution  as  the  flow  on  the 
artificial  arc  is  reduced  to  zero. 


6.3  Parametric  Node  Incremental  Flow  Algorithm 

l 

Given  a  network  with  arcs  having  capacity,  cost  and  gain 
parameters,  with  source  and  destination  nodes  defined,  and  required 
quantities  of  flow  specified  for  each  source  node  and  destination 
node,  the  problem  is  to  determine  the  arc  flows  with  the  minimum 
total  cost.  The  flows  must  satisfy  the  conservation  of  flow  except 
at  the  super-source.  While  the  flows  into  each  destination  must 
equal  the  required  values,  the  super-source  may  supply  any  amount  of 
flow.  The  problem  may  be  stated 

m 

Min  E  h.  f,  (6.8) 

k-1  K  * 

S.T.  -  E  a.  f.  +  E  f.  -  b.  for  ie  N-s  (6.9) 

^“ti  ke  M0i 

°  <  fk  <  ck  •  (6.10) 

The  parametric  node  incremental  algorithm  is  a  dual-node-infeasible 
approach  applied  to  a  specialized  network.  The  specialized  network 
has  a  single  source,  but  perhaps  multiple  destination  nodes 


I. 


Min 


m 

s  Vk 


k-1 


S.T.  -  Z 


k£ 


Vk 


for  i £  N-s 


k£  M, 
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where  b  b^  for  1  eN  -  1.  An  intermediate  optimum  is  a  flow 

and  node  potentials  t that  are  an  optimum  solution  to  the  above 
primal  problem  and  its  associated  dual. 

The  general  approach  begins  with  an  intermediate  optimum  Fn 


for  some  value  of  output  flow  at  each  of  the  sinks.  Next  a  shortest 


path  network  is  constructed,  D  ,  that  determines  for  each  node  the 

s 


minimum  cost  per  unit  of  additional  flow  to  the  sink  nodes  and  the 


path  over  which  that  flow  may  be  obtained.  The  output  flow  at  each 


0 


destination  node  is  increased  parametrically  in  Dg  along  the 


minimum  cost  trails  defined  to  all  the  destinations  until  an  arc 


becomes  capacitated.  This  flow  augments  Fq  to  Fj.  Fj  is  also  an 


Intermediate  optimum.  A  new  augmenting  network  is  constructed,  D 


and  the  output  flow  is  again  augmented  to  obtain  F£.  This  process 
continues  iteratively  until  Che  desired  output  flows  at  the 


destination  nodes  are  met 


Throughout  the  algorithm  because  the  values  of  tt  ^  are 
determined  by  a  shortest  path  algorithm,  the  solution  will  remain 
dual  feasible.  The  solution  is  always  basic  as  at  most  n-1  arcs 
will  have  flows  strictly  between  their  bounds,  and  these  arcs  will 
form  a  basis  network.  The  solution  will  satisfy  conservation  of 
flow  and  complementary  slackness.  The  one  requirement  that  is 
violated  until  the  final  iteration  is  the  requirement  for  the  flow 
at  each  of  the  destination  nodes  to  be  satisfied.  These  flows  start 
at  zero  and  iteratively  Increase  until  they  reach  their  required 
values.  At  this  point  the  algorithm  stops  with  the  optimum 
solution.  The  algorithm  has  the  following  basic  steps: 


1.  Find  an  initial  least-cost  spanning  tree. 

2.  If  all  sink  external  flows  are  met,  stop 
with  the  optimal  solution. 

3.  For  each  sink  node,  find  the  flow 
augmenting  trail.  Parametrically  increase 
the  flow  on  all  trails.  Determine  the 
maximum  allowable  flow  increment  and 
determine  the  arc  to  leave  the  network. 

4.  Find  an  arc  to  enter  the  network.  If  none 
exists,  stop .  There  is  no  feasible 
solution. 

5.  Change  the  basis  by  deleting  the  leaving 
arc  and  inserting  the  entering  arc.  Modify 
the  dual  variables  to  satisfy  complementary 
slackness.  Return  to  step  2. 
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Again  before  stating  the  algorithm  in  complete  detail,  each  step  of 
the  outline  will  be  discussed. 

6.3.1  Obtaining  the  Initial  Basis  Tree 

In  order  to  apply  the  concepts  developed  above  and  the  ones 
presented  in  Chapter  5,  we  need  to  modify  the  network  as  presented. 
Consider  the  sample  network  shown  in  Figure  6.5.  We  create  a 
super-source  and  arcs  from  the  super-source  to  each  source  node  with 
parameters: 


\  "  bsi 
\  *  1  * 


This  transformation  is  accomplished  by  the  subroutines  READNDE  and 
ORIGSG.  The  revised  network  is  presented  in  Figure  6.6.  Once  the 
transformed  network  is  obtained,  we  use  a  generalized  Dijstra's 
algorithm  as  documented  in  [46, page  283]  to  obtain  in  initial  basic 
solution.  We  then  use  the  subroutine  TREINT  [46, page  121]  to 
initialize  the  parameters  of  the  basis  tree.  Since  a  generalized 
network  problem  may  have  arc  gains  greater  than  one  or  have  negative 


The  Augmented  Network 
Figure  6.6 
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arc  costs,  Dijkstra's  algorithm  may  not  provide  an  optimal  shortest 
path  tree.  Therefore,  we  now  use  the  algorithm  PSHRTG  [46, page  285] 
to  find  the  initial  least-cost  spanning  tree.  At  this  point  we  have 
an  initial  optimal  (all  d^  0) ,  basic,  node  infeasible  solution  to 
the  original  problem. 

6.3.2  Determining  the  Arc  to  Leave  the  Basis 

As  in  PSENG,  we  now  use  the  PATH  subroutine  to  determine  the 
flow  augmenting  paths.  In  this  case,  however,  the  candidate  list 
contains  node  numbers,  and  the  node  parametric  parameters  are  the 
negative  of  the  respective  destination  node  external  flows.  For 
each  destination  node  we  find  the  flow  augmenting  trail.  This  will 
create  a  list  of  arcs  (LISA),  a  list  of  nodes  (LISN),  and  a  gamma 
list  (G).  Using  these  lists,  the  subroutine  MFLOWG  [46, page  266] 
will  determine  the  maximum  flow  change  or  reference  number  on  these 
trails  and  the  limiting  arc  k^. 

At  each  iteration  we  keep  a  record  of  the  cumulative 
reference  number.  This  is  merely  the  sum  of  the  reference  numbers 
to  this  iteration.  It  is  possible  that  the  reference  number  found 
for  a  given  iteration,  1,  when  added  to  the  cumulative  reference 
number  from  the  previous  iteration  1-1,  will  be  greater  than  1. 

This  means  a  reference  number  of  (i  -  cumulative  reference  number) 


will  cause  Che  unsatisfied  external  flows  Co  be  satisfied  exactly. 

At  this  point  we  would  augment  the  flows  and  terminate  with  the 
optimal  solution. 

6.3.3  Determining  the  Arc  to  Enter  the  Basis 

We  use  the  subroutine  FINDNDE  to  determine  the  entering  arc 
k^.  If  no  arc  can  be  found  to  enter  the  basis,  no  feasible  solution 
exists  to  the  proposed  problem. 

6.3.4  Changing  the  Basis 

Updating  the  basis  tree  representation  and  the  dual 
variables  Is  accomplished  by  the  subroutine  PIV0T1G. 


6.3.5  The  Complete  Algorithm 


Use  a  variant  of  Dijstra's  shortest  path 
algorithm  to  find  the  initial  least-cost 
spanning  tree.  Use  the  arcs  in  the 
spanning  tree  to  define  a  set  of  node 
labels  it  ^  such  that  ir  g  ■  0  and  tt  j  »  (  ir  ^ 
+  1^)  /  afc  for  all  arcs  k(i,j)e  Hj.. 

Use  a  triple  labeling  scheme  to  label  and 
store  the  least  cost  spanning  tree. 

If  the  cumulative  reference  number  (CUMRNO) 
equals  1.0,  stop  with  the  optimal  solution. 

For  each  demand  node  find  the  flow 
augmenting  trail.  Determine  the  maximum 
reference  number,  MF,  on  these  trails.  If 
the  maximum  MF  is  such  that  the  CUMRNO  +  MF 
>  1,  set  MF  -  1  -  CUMRNO  and  terminate  the 
algorithm  after  this  step.  Determine  the 
arc  on  these  trails  which  is  at  capacity 
and  designate  it  the  leaving  arc  k^. 

Augment  the  flows. 


6 


Delete  the  arc  k^  from  the  basis  tree, 
partitioning  the  network  into  two  subtrees 
defined  by  the  set  N^  and  the  set  N^. 

Search  for  an  admissible  arc  that 
originates  or  terminates  in  N^  and  connects 
to  the  set  or  an  admissible  arc  that 
originates  and  terminates  in  N£.  If  no 
such  arc  exists,  stop;  there  is  no  feasible 
solution  to  the  problem.  Otherwise  select 
the  arc  with  the  minimum  and  designate 
it  as  the  entering  arc  kg. 

7.  Perform  a  basis  change  which  determines  the 
new  least'cost  spanning  tree.  The  basis 
change  is  accomplished  by  deleting  k^  and 
adding  kg.  Update  the  dual  variables  and 
return  to  step  4. 


The  algorithm  PARANDE  implements  these  steps.  The  flow  charts  for 
the  FORTRAN  coding  of  the  algorithm  can  be  found  in  Appendix  2. 

Figure  6.7  presents  an  example  application  of  the  algorithm 
PARANDE.  Figure  6.7(a)  shows  the  original  problem  data.  In  each  of 


the  Figures  6.7(b)  through  6.7(h),  the  entering  arc  is  indicated  by 
a  dashed-dot  line,  the  leaving  arc  by  an  X,  the  basic  arcs  by  solid 
lines  and  the  nonbasic  arcs  at  capacity  by  dashed  lines.  The 
unsatisfied  external  flows  are  also  shown  at  each  step.  Figures 
6.7(b)  through  6.7(g)  show  the  intermediate  steps  in  the  solution 
process.  In  Figure  6.7(g)  no  leaving  arc  is  shown,  as  the  maximum 
flow  augmentation  found  is  greater  than  the  maximum  flow 
augmentation  necessary  to  satisfy  the  external  flows.  The  algorithm 
terminates  at  this  point  and  the  optimal  solution  is  shown  in  Figure 
6.7(h). 

6.4  Computational  Analysis  of  PARARC  and  PARANDE 

To  demonstrate  the  feasibility  of  using  parametric  analysis 
for  Incremental  flow  algorithms,  PARARC  and  PARANDE  were  tested 
against  their  immediate  predecessor  INCREMG  [46, page  314].  INCREMG, 
in  contrast  to  the  parametric  approach,  only  augments  flow  along  a 
single  flow  augmenting  path.  In  the  parametric  approach  flow  is 
increased  in  multiple  flow  augmenting  paths.  The  computations  were 
conducted  to  determine  if  these  new  algorithms  would  be  more 
efficient  than  INCREMG.  The  computational  times  were  also  compared 
against  two  primal  codes,  NETG  (38)  and  PGAINS  (46).  We  should  note 
that  both  PGAINS  and  INCREMG  are  modular  and  pedagogical  in  design. 
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Example  Application  of  PARANDE  Algorithm 
Figure  6.7(cont) 
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Figure  6.7(cont) 
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Example  Application  of  PARANDE  Algorithm 
Figure  6.7(cont) 
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and  this  may  influence  the  computational  times  recoreded. 

A  set  of  20  test  problems  (see  Table  6.1)  included  in  the 
computational  analysis  were  generated  by  a  random-network  generator 
NETGNG.  While  NETGNG  has  no  specific  reference,  it  is  in  use  by  the 
Center  for  Cybernetic  Studies  at  The  University  of  Texas  at  Austin. 
NETGNG' s  immediate  predecessor  for  pure  networks  is  NETGEN  (53). 
NETGEN,  as  part  of  its  procedure,  builds  a  skeleton  network 
connecting  sources  to  sinks.  Capacities  are  assigned  to  all 
skeleton  branches  to  assure  the  existance  of  a  feasible  solution 
which  uses  all  supplies  and  meets  all  demands.  The  sum  of  supplies 
equals  the  demand.  Nonskeleton  branches  are  then  randomly 
generated.  The  primary  modification  to  this  procedure,  used  by 
NETGNG,  for  networks  with  gains  is  also  to  assign  nonskeleton 
branches  a  gain  parameter  randomly  selected  from  a  specific  range. 
The  skeleton  branches  are  given  unity  gains.  This  guarantees  the 
existance  of  a  feasible  solution. 

The  problem  set  includes  10  transportation  and  10 
transshipment  networks.  The  number  of  nodes  varies  from  60  to  100 
and  the  number  of  arcs  from  100  to  600.  The  cost  on  the  arcs  ranges 
between  1  and  100.  All  problems  have  100%  capacitated  arcs,  with 
minimum  lower  bound  of  5  and  maximum  upper  bound  of  10.  Arc  gains 
range  between  0.5  and  1.5.  All  problems  were  generated  by  using  a 
seed  of  12345678  and  80%  of  the  arcs  were  assigned  the  highest 


Problem  Number  Number  Number  Number  Cost  Total  Capacity  Bound  Gains 
Number  of  of  of  of  Range  Supply  Min  Max  Min 

Nodes  Sources  Sinks  Arcs 


Table  6.1 


possible  cost. 

The  five  algorithms  were  run  on  a  Cyber  170/750A  computer  at 
The  University  of  Texas  at  Austin.  All  codes  were  tested  on  FTN4 
(OPT  »  2)  compiler.  Tables  6.2  and  6.3  display  the  results  from 
NETG,  P GAINS,  INCREMG,  PARARC,  and  PARANDE.  All  solution  times 
include  only  the  time  to  optimize  the  problem  excluding  times  to 
input  the  data  and  output  the  results. 

It  can  be  concluded  based  on  the  limited  sample  results 
shown  in  Tables  6.2  and  6.3  that  both  primal  codes  completely 
dominated  the  dual-incremental  codes  by  between  two  to  five  times. 
However,  PARANDE  is  more  efficient  than  either  PARARC  or  INCREMG. 

For  the  10  transportation  problems  PARANDE  is  1.41  times  faster  than 
PARARC  and  1.67  times  faster  than  INCREMG.  For  the  10  transshipment 
problems  PARANDE  is  1.27  times  faster  than  PARARC  and  1.67  times 
faster  than  INCREMG.  Overall  PARANDE  is  1.31  times  faster  than 
PARARC  and  1.67  times  faster  than  INCREMG. 

It  can  also  be  concluded  that  PARARC  is  more  efficient  than 
INCREMG.  For  the  transportation  problems  PARARC  is  1.18  times 
faster  than  INCREMG  and  for  the  transshipment  problems  PARARC  is 
1.31  times  faster  than  INCREMG.  Overall  PARARC  is  1.27  times  faster 
than  INCREMG. 

The  major  objective  of  this  computational  study  was  to 
demonstrate  the  feasibility  of  applying  parametric  analysis  to 


Problem NETG PGAINS  INCREM6  PARARC  PARANDE  Objective 
Number  I  Function 


0.316  0.968 


0.960  1.425 


0.942  1.587 


0.657  1.761 

1.043  1.811 


0.632  2.029 


1.260  2.594 


0.131  0.564 


0.565  1.138 


1.032 


1.966 


3.863 


4.155 

5.157 


6.183 


10.682 


0.670 


3.881 


0.958 


2.087 


3.965 


8.066 


0.687 


3.151 


0.622 


1.383 


2.697 


3.931  2.635 
4.717  4.125 


5.542  3.745 


5.487 


0.487 


2.289 


844765.0 


743281.6 


801090.7 


840359.3 
746985.6 

812475.3 
757218.0 


831556.8 


715154.7 


0.693  1.507 


4.775 


2.689 


1.866 


743944.5 


Solution  Times  for  the  Transportation  Problems 
Table  6.2 


PRIMAL  METHODS  I 

DUAL  METHODS  1 

Problem 

Number 

||||||^ 

INCREMG 

PARARC 

PARANDE 

Obj ective 
Function 

11 

0.453 

1.146 

4.738 

4.664 

3.666 

2357360.3 

12 

0.828 

1.449 

7.096 

5.408 

4.276 

2085017.7 

13 

0.936 

1.812 

10.432 

7.243 

5.662 

2203213.5 

14 

1.516 

2.310 

18.143 

13. 394 

10.959 

2082679.7 

15 

0.915 

1.296 

5.068 

5.471 

4.223 

2230957.7 

16 

0.913 

1.587 

8.100 

7.072 

5.655 

2473342.4 

17 

1.143 

2.929 

15.195 

9.785 

7.746 

2117899.7 

18 

0 . 969 

1.886 

8.440 

6.375 

4.696 

2055368.9 

19 

0.986 

2.139 

12.480 

10.116 

7.698 

1951960.6 

20 

1.616 

2.925 

21.410 

15.285 

12.053 

2169286.6 

Solution  Times  for  the  Transshipment  Problems 
Table  6.3 


dual-incremental  methods.  While  the  codes  tested  were  not  as 


efficient  as  primal  codes,  they  were  considerably  more  efficient 
than  the  previous  generation  of  dual-incremental  codes.  This  does 
indicate  an  area  of  promising  research  which  will  be  discussed  in 
Chapter  7. 

6.5  Other  Applications  of  Parametric  Analysis 

In  this  section  we  will  describe  four  additional  problem 
classes  where  parametric  sensitivity  analysis  for  generalized 
networks  might  be  applied,  with  what  we  feel  would  be  with  a  great 
deal  of  success. 

The  first  is  the  cash  management  problem  (9)  which  is 
concerned  with  optimally  financing  net  cash  outflows  and  investing 
net  cash  inflows  of  a  firm  while  simultaneously  determining  payment 
schedules  for  incurred  liabilities.  This  could  be  formulated  as  a 
multi-period  transshipment  model  with  m  "warehouses"  corresponding 
to  different  sources  of  cash  (cash  sales,  short-term  securities, 
lines  of  credit,  etc.)  and  the  n  "markets"  corresponding  to 
different  uses  of  cash  (payments  on  accounts,  notes  payable,  etc.). 
Cash  could  be  both  a  source  and  a  use,  thus  allowing  transshipments 
across  periods.  Yields  on  securities,  interest  rates  on  loans,  and 
discounts  constitute  unit  costs  of  the  problem.  While  the  model 


Cakes  into  account  payment  schedules,  financing  and  security 
transactions,  it  Ignores  the  Important  problem  of  determining  the 
minimum  balance.  Cash  deposits  in  excess  of  some  absolute  minimum 
cash  balance  requirement  prove  valuable  since  they  improve  the 
firm's  credit  rating  and  the  bank's  goodwill.  But  deposits  also 
mean  lost  revenue  from  alternative  investments  in  securities.  Thus 
the  manager  would  like  to  know  the  effect  on  the  optimum  total  cost 
from  maintaining  different  levels  of  cash  balance,  so  as  to  arrive 
at  a  subjective  decision  regarding  its  magnitude.  This  can  be 
easily  accomplished  by  parametric  analysis  of  the  cash  balance 
(external  flows). 

A  second  major  area  is  the  multi-period  multi-region 
capacity  expansion  problem.  Balachandran,  et  al.  (9),  describe  this 
problem.  A  product  can  be  manufactured  in  m  regions  and  is  required 
by  n  markets  for  each  of  the  periods  K  *  {1,2,...,T}.  Ue  know  that 
the  Initial  demand  and  the  demand  in  a  market  can  be  satisfied  by 
production  and  shipment  from  any  region  but  must  be  met  exactly 
during  each  time  period  (no  backlogging  or  inventorying) . 
Balachandran  notes  that  each  period  is  a  transportation  problem,  and 
the  problem  for  period  t  differs  from  the  problem  in  period  t-1  by 
only  the  right-hand  side  constraints.  Thus  we  can  iteratively  solve 
the  problem.  Starting  in  period  1  and  using  the  solution  from 
period  1,  we  use  parametric  analysis  to  get  the  solution  for  period 
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2.  We  can  continue  In  a  like  manner  to  get  the  solution  for  period 
T. 

A  third  major  area  would  be  In  the  area  of  branch-and-bound 
algorithms.  It  often  happens  that  the  subproblems  are  generalized 
network  flow  problems.  Consequently,  as  we  move  along  the 
branch-and-bound  tree,  we  may  obtain  the  optimal  solution  to  the 
(k+l)st  subproblem  by  post-optimization  of  the  kth  subproblem, 
rather  than  Inefficiently  having  to  solve  the  (k+l)st  problem  from 
the  start.  Some  of  the  potential  areas  within  branch-and-bound 
might  deal  with  optimal  facility  location,  assignment  of  sources  to 
users,  lock-box  decision  models,  or  the  one-machine  job  scheduling 
problem. 

A  final  area  where  this  method  might  make  solution  methods 
more  efficient  is  in  the  area  of  stochastic  network  flow  problems 
with  recourse  (45).  In  this  problem  one  really  has  a  master  problem 
and  a  series  of  subproblems.  Each  of  the  subproblems  is  a  network 
where  the  right-hand  side  is  being  varied  at  each  iteration.  It 
seems  that  the  parametric  sensitivity  algorithm  would  be  of  value  to 
solve  these  types  of  problems,  as  the  solution  for  iteration  k  could 
be  used  to  start  iteration  k+1  by  just  varying  the  capacity 
constraints. 

While  the  four  example  application  areas  do  not  exhaust  the 
potential  areas  of  applications  for  the  parametric  analysis  method 


for  generalized  networks,  we  hope  that  these  examples  may  spur 
further  applications  along  these  lines. 


CHAPTER  7 


CONCLUSIONS  AND  RECOMMENDATIONS 


7.1  Summary 

In  a  deterministic  generalized  network  after  the  Initial 
model  Is  constructed,  the  operations  research  scientist  must  still 
gather  the  data  to  completely  fill  in  the  parameters  of  the  model. 
The  uncertainty  and  Inaccuracy  of  the  Initial  data  collection  is 
well  documented.  When  the  operations  research  scientist  has 
presented  these  results  to  management,  the  scientist  has  many  times 
met  with  failure  because  management  challenges  the  original  problem 
data.  The  real  strength  of  linear  programming  has  been  that  the 
operations  research  scientist  has  been  able  to  use  sensitivity 
analysis  or  parametric  programming  to  perturb  the  Initial  problem 
data  to  answer  management's  "what  if"  questions. 

A  logical  extension  of  parametric  programming  of  linear 
programs  would  be  to  apply  these  well  known  concepts  to  a 
generalized  network  format.  The  goal  of  this  research  was  to 
develop  the  individual  algorithms  to  accomplish  this  both  by 
parametric  analysis  for  a  capacity  vector,  as  well  as  parametric 
analysis  of  the  right-hand  side  or  external  flow  vector.  In  order 


to  accomplish  this  goal,  an  algorithm  to  determine  the  flow 
augmenting  trails  (TRAIL2) ,  an  algorithm  to  determine  the  maximum 
flow  augmentation  possible  (MFL0G2),  an  algorithm  to  determine  the 
entering  arc  (FIND2G),  and  finally  an  algorithm  to  change  the  basis 
(PIV0T1G)  were  developed. 

A  computer  package  was  developed  which  integrated  these 
algorithms  with  a  report  generator.  The  report  allows  the  manager 
to  see  the  optimal  flows  at  each  iteration  and  examine  how  the  total 
cost  varies  with  each  basis  change. 

•  An  outgrowth  of  this  research  has  been  the  development  of 
two  new  variants  of  a  dual-incremental  flow  approach  to  the 
generalized  network  flow  problem. 

7.2  Contributions 

As  described  in  Chapter  2,  much  work  historically  has  been 
spent  on  parametric  programming  for  linear  programs.  This  is  in 
contrast  to  a  lack  of  significant  work  in  the  area  of  parametric 
programming  for  the  generalized  network  flow  problem.  There  are  two 
major  contributions  of  this  research.  First,  the  methods  developed 
in  this  study  have  not  been  seen  before  in  the  literature  and  they 
are  a  significant  tool  for  the  practical  operations  research 
analyst.  The  methods  allow  the  analyst  to  perturb  an  initial 
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uncertain  data  set  to  analyze  how  these  perturbations  will  effect 
recommendations  to  management.  The  computer  code  developed  allows 
the  analyst  to  present  a  graphical  report  on  how  the  objective 
function  varies  as  the  parameters  are  varied.  This  leads  to  an  easy 
visualization  of  how  much  to  vary  the  parameters  or  how  much  the 
parameters  may  be  varied  before  an  adverse  effect  is  likely  to 
happen  to  the  company. 

While  the  first  contribution  is  of  great  practical  value, 
the  second  outgrowth  of  this  research  has  been  the  development  of 
two  new  dual-incremental  codes.  While  these  codes  are  not 
numerically  superior  to  primal  methods,  they  do  extend  the  theory  of 
parametric  analysis  to  solve  generalized  network  flow  problems. 
Before  this  research,  parametric  analysis  as  related  to  operator 
theory  was  used  only  to  solve  the  generalized  transportation  problem 
(5, 6, 7, 8). 

7.3  Recommendations  for  Further  Research 

There  are  three  areas  which  could  lead  to  further  research. 
The  first  logical  extention  would  be  to  develop  parametric 
programming  codes  for  changes  in  arc  costs  or  gain  parameters.  This 
would  allow  the  operations  research  analyst  to  vary  all  the 
parameters  of  the  network  to  do  a  complete  analysis  of  the  proposed 


model. 

The  second  area  would  be  to  continue  further  testing  to  make 
PARARC  and  PARANDE  computationally  more  efficient.  In  (67)  Schmitt 
developed  a  pure  dual-incremental  code  which  is  relatively 
competitive  with  known  primal  codes.  This  leads  to  speculation  that 
dual-incremental  codes  for  generalized  networks  might  also  be 
competitive  if  they  were  modified  appropriately.  The  principle  time 
consumer  for  a  dual  method  is  the  search  for  the  arc  to  enter  the 
basis.  This  usually  requires  a  search  in  each  iteration  proportional 
to  the  number  of  arcs.  In  a  private  communication  (58)  Matsumoto 
suggests  that  one  only  need  search  through  those  admissible  arcs 
which  originate  or  terminate  in  the  Nj  set.  Another  potential  way 
to  reduce  computational  time  is  to  implement  PARARC  and  PARANDE  with 
the  augmented  threaded  index  method.  Glover,  fUlngman,  and  Stutz 
report  a  10Z  reduction  in  computation  time  using  preorder  traversal 
lists  for  the  basis  representation  rather  than  the  triple  label 
method. 

A  third  way  to  improve  PARARC  and  PARANDE  would  be  the  way 
the  initial  basis  tree  is  found.  While  the  basis  initialization 
only  consumes  2-3Z  of  the  total  optimization  time,  it  is  possible 
that  the  algorithm  THRESH  presented  by  Glover,  et  al.  (39)  may  be 
able  to  reduce  this  even  further.  A  fourth  possible  way  to  save 
time  is  to  eliminate  the  use  of  mirror  arcs  in  a  way  similar  to 
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Schmitt's  work  (67).  The  comparison  to  check  if  arc  k  is  greater 
than  zero  is  done  numerous  times  throughout  the  complete  algorithm. 
Without  counting  the  number  of  times  this  logical  check  is  done 
after  the  computations  with  do  loops  are  applied,  the  author  has 
roughly  counted  over  40  individual  logical  comparisons  of  k  >  0.  A 
fifth  possible  way  to  lower  the  optimization  time  is  to  continue  to 
Increment  the  flow  on  those  arcs  on  the  tree  until  all  paths  from 
the  source  to  the  destinations  are  capacitated.  Then  by  using  a 
modified  PSHRTG  a  new  basis  tree  could  be  found.  This  could 
potentially  cut  down  on  the  number  of  dual  Iterations,  thus  greatly 
reducing  computational  time. 

A  final  area  of  research  Interest  would  be  to  apply  the 
parametric  methods  to  multicommodity  flow  problems.  There  has  been 
much  recent  Interest  in  this  area  (16,34,35,51,60,65).  It  seems 
possible  that  a  two  commodity  flow  problem  in  the  form: 


min  H  F  +  H  F 
y  y  xx 

S.T.  AF  -  b 

y  y 

AF  -  b 
x  x 

0  <  F  <  GF 
x  —  y 

Fy  2  0 


could  be  solved  by  dividing  this  problem  into  a  master  problem  and  a 


I 


L 


subproblem.  In  Che  subproblem  Che  concepCs  of  parameCric  analysis 
could  be  used  Co  modify  Che  rlghc-hand  sides  and  Chen  pass  Che 
required  daca  back  Co  Che  mascer  problem. 


MAINARC 


PURPOSE:  To  read 
data,  initialize 
the  basis,  and 
call  parametric 
are  incremental 
algorithm  to 
generalized  net- 
work  problem. 


1. (DATA)  Read 
problem  data  from 
input  source. 

2.  (INITIAL)  Use 
DSHRTG  to  find  the 
shortest  path  from 
slack  node  to 
super-sink  node. 

If  no  path  exists 
(NP  *1),  stop. 
Otherwise  use 
TREINT  to  con¬ 
struct  pointer  rep¬ 
resentation  of 
tree.  Call  PSHRTG 
to  find  shortest 
path  for  general¬ 
ized  net  work.  If 
UNB  “1,  an  unbound¬ 
ed  solution  exists. 
Otherwise  go  to  3. 


c 


DATA 

(  READGl  ^ 

INITIAL 

(  D5HRTG(SN,TN,NP)  ) 

NP  -1 

STOP 

(  TREINT(N)  ) 

(  PSHRTG(UNL.J)  ) 

Vy  UNB  -1 

STOP 

- > ALGORITHM 

3. (ALGORITHM)  Call 
parare  to  solve 
problem.  If  ST  -  1 
there  is  no  feasi¬ 
ble  solution.  GEN- 
OUT  prints  solution 
results. 


Main  Program  MAIN ARC 


PURPOSE:  To  accept 
arc  data  itea  and 
•tore  It  in  an  arc 
list  ordered  by  as¬ 
cending  origin 
node. 

1. ( INITIAL)  If 
first  call  to 
0RIGSG1,  set  all 
node  pointers  to  1. 
Otherwise  go  to  2. 

2.  (MOVE)  Increase 
M  by  one.  Increase 
all  node  pointers 
greater  than  I  by 
one.  Move  all  arcs 
above  the  new  en¬ 
try  one  index  high¬ 
er  on  list. 

3.  (ARC)  Insert  arc 
in  last  position 
alloted  to  node  I. 
Modify  the  upper 
bounds  of  the  arcs 
and  fixed  external 
flows  to  account 
for  the  arc  lower 
bound • 


Subroutine  0RIGSG1 


a 


PURPOSE:  To  read 
and  score  node  and 
arc  data  for  che 
minimum  cose  flov 
problem. 

1.  (INITIAL)  Init¬ 
ialise  che  number 
of  arcs  to  zero. 
Create  super-slack 
node  and  super- 
sink  node  (xs). 

Sec  all  external 
flows  to  zero. 

2.  (NODE)  Read  a 
node  data  item. If 
fixed  external  flow 
is  positive, put 
fixed  flow  in  stor¬ 
age.  Otherwise,  add 
up  the  cumulative 
amount  of  demands 
and  create  arc  from 
demand  node  to  xs. 

If  slack  external 
flow  is  zero  go  to 

2.  Otherwise  create 
a  slack  arc  and 
store  data  in  cor¬ 
rect  position  in  arc 
lists.  After  all 
node  data  has  been 
read  (1-0).  Set  sup¬ 
er-sink  external 
flow  to  cumulative 
demands.  Go  to  3. 

3.  (ARC)  Read  an  arc 
data  item.  If  1*0, 


READG1 


INITIAL 


READ, N,M : -0 .SLACK : -N+2 , XS : N+l , XM : -0, XX: -0 
FOR  I:-l  TO  N 
//  B(I):-0 


READ  I,BF,BS,CS 


I  -  0 


fa 


B(I):-BP  XX : -XX+BF , J : -XS , LOWER : -0 , 

UPPER :-0, COST : -BIG+1 , GAIN : -1 

f 0RICSC1 ( I , J ,L0UER, UPPER ,  \ 
V. _ COST,  GAIN) _ J 

K  BS  -  0  / 


Y\^ _ BS  >  0 _ / 

J: -I,  I -.-SLACK,  J: -SLACK,  LOWER: -C 

LOWER  :-0,  UPPER  :-BS  UPPER:— BS, 

COST: -CS, GAIN: -l  COST :-CS, GAIN :-l 


LI 


0RIGSG1(I,J, LOWER, UPPER, COST.GAIN) 


- >N0DE 


B(XS):-XX,BF:-0 


- >ARC 


Subroutine  READG1 
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PURPOSE:  To  im- 
plement  Che  para- 
mecric  arc  flow 
augmentation  ap¬ 
proach  for  a  gener¬ 
alized  network. 

1.  (RESET)  Set  the 
costs  on  all  arcs 
terminating  at  the 
super-sink  node  to 
zero,  except  the 
artificial  arc. 

2.  (INITIAL)  Init¬ 
ialize  all  para¬ 
meters  to  zero. 

3.  (OPTIMAL)  If  the 
flow  on  the  artifi¬ 
cial  arc  equals  0, 
return  with  opti¬ 
mal  solution. 
Otherwise  go  to  4. 

4.  (CHECKIN)  Init¬ 
ialize  the  T's  to 
zero.  Initialize 
arc, node,  and  cycle 
counters  to  zero. 

5.  (LEAVE)  Por  each- 
arc  on  the  candi¬ 
date  arc  list ,  de¬ 
termine  if  the  arc 
is  nonbaslc.  If  it 
is  and  its  para¬ 
metric  parameter  is 
less  than  zero, 
calculate  a  refer- 


CHECKIN 


1/ 

II 

FOR  I:-l  TO  N 

G(I):-0,IK(I):-0 

IC:-0,CK:-0,Z:-0,ME:-99999,NK:-0 

- > LEAVE 

Subroutine  PARARC 


ence  number.  Score 
Che  smallest  refer¬ 
ence  number  ME.  If 
Che  paramecer  la 
negative,  is  sec  IE 
co  che  cerminal 
node  of  XX,  other¬ 
wise  IE  Is  che 
origin  node  of  XX. 
Call  PATH  which 
creates  LISA.LISN 
and  RC  Uses  for 
all  flow  augment¬ 
ing  trails.  If  a 
cycle  has  been 
found  in  PATH,  call 
CYCLEG  co  updaCe  1 
values  for  all 
nodes  on  Che  cycle. 
Call  MFL0G2  Co  de¬ 
termine  Che  maxi¬ 
mum  reference  num¬ 
ber  and  arc  which 
reaches  capacity. 

If  XL  •  0  any  fur¬ 
ther  iterations 
trill  drive  a  cap- 
city  below  zero.  If 
che  reference  num¬ 
ber  is  zero,  this 
is  a  degenerate 
iteration.  Other¬ 
wise  update  flows 
on  basic  arcs  and 
flows  and  capaci¬ 
ties  on  candidate 


LEAVE 


I:»l  TO  XM 


XXX: -XX( I) , ZI : O (XXX) , ZJ : -T(XXX) 

Vv  PB(ZJ)  -  XXX  OR  PB(ZI)  -  -XXX 

X.  P(XXX)  <  0 

Ml : -F(XXX) /-PARA(I)  XXE:-XXX 
PX:— PARA(I)  PX:-PARA(I) 

XXE:— XXX  IE:  -O(XXX) 

IEs-T(XXX)  JE: -T(XXX) 

JE:-0(XKK) 

X  Ml^ME  /N  [ 

ME: -Ml, NX: -XXX  I  “ 

I J :  -IE ,  TM :  -1 ,  MF :  -ME ,  ITER :  -ITER+1 

C  PATH( I J , TM , XXE , XXX, PX, IX, CX) 

IJ:-JE,TM:—  1 

(  PATH(IJ,TM.XXE,XRX,PX,IX, CX) 


Z  -  0 


MFL0G2 


XL  -  0 


CYCLEG 


PL0WG2 


ITERD: - 
ITERD+1 

- >ENTER 


FL0WG2 


Subroutine  PARARC(cont) 


ENTER 


6.  (ENTER)  Call 
FIND ARC  (same  as 
FIND2G  except  it 
has  been  modi¬ 
fied  for  FOR¬ 
TRAN  coding)  which 
will  determine  KE. 
If  KE  -  0  there  is 
no  feasible  solu¬ 
tion  to  the  prob¬ 
lem  and  return* 
Otherwise  go  to  7 . 

7.  (CHANGE)  Change 
the  basis  and  up¬ 
date  dual  vari¬ 
ables*  Go  to  3. 


Subroutine  PARARC(cont) 


FL0MG2(IC,MF) 


PURPOSE:  To  change 
flow  on  each  basic 
arc  by  MF  and  Co 
update  candidate 
arc  flows  and 
capacities. 


INITIAL 


r\^  ic  -  o 

// 

If 

II 

n 

n 

// 

u 

// 

// 

// 

FOR  L:-l  TO  IC 

K:-LISA(L),J:-LISN(L) 

K  -  0 

K _ v>° _ . 

F(K):- 

F(K)-WF*G(J)/ 

A(K) 

F(-K):- 

F(-K)-MF*G(J) 

- >PARAMETRIC 

PARAMETRIC 


FOR  JJJ:-1  TO  XM 

- 

ft  K:-XK(JJJ),C(K):-C(K)*MF*PARA(JJJ), 
U  F(K):-F(K)+MF*PARA(JJJ) 


ETURN 


Subroutine  FL0WG2 


PURPOSE:  To  read 
data,  initialize  - 
the  basis,  and 
call  parametric 
node  incremental 
algorithm  to 
generalized  net¬ 
work  problem. 


1. (DATA)  Read 
problem  data  from 
input  source. 

2.  (INITIAL)  Use 
DSHRTG  to  find  the 
shortest  path  from 
slack  node  to  all 
destinations. 

If  no  path  exists 
(NP  -1),  stop. 
Otherwise  use 
TREINT  to  con¬ 
struct  pointer  rep¬ 
resentation  of 
tree.  Call  PSHRTG 
to  find  shortest 
path  for  general¬ 
ized  net  work.  If 
UNB  «l  an  unbound¬ 
ed  solution  exists. 
Otherwise  go  to  3. 


(  MAINNDE  ) 


DATA 

(  READNDE  ) 

INITIAL 

(  DSHRTG(SN.TN.NP)  ) 
- ' 

Y\  NP  -1 

(  TREINT(N)  ) 

(  PSHRTG (UNB , J )  ) 

UNB  -1 

[stop  STOP 

- > ALGORITHM 

ALGORITHM 

(  PARANDE  ) 

K.  ST  •  1 

\ _ _ _ I _ _ _ Z 

SO  FEASIBLE  SOLUTION 

(  GENOUT  ) 

3. (ALGORITHM)  Call 
parande  to  solve 
problem.  If  ST  •  1 
there  is  no  feasi¬ 
ble  solution.  GEN- 
OUT  prints  solution 
results. 


Main  Program  MAINNDE 


INITIAL 


READNDE 


PURPOSE:  To  read 
and  store  node  and 
arc  data  for  Che 
ainiaua  coat  flow 
problea . 

1.  (INITIAL)  Init¬ 
ialize  the  nuaber 
of  arcs  to  zero, 
total  deaand  to 
zero,and  nuaber  of 
paraetric  nodes  to 
zero.  Set  SLACK  - 
nodes  +  1. 

2.  (NODE)  Read  a 
node  data  itea.If 
the  itea  is  blank 
go  to  3.  If  the 

f lxed  external  flow 
is  positive  put  it 
in  storage  .  Other¬ 
wise  set  up  para- 
aetric  node  list  and 
keep  track  of  total 
deaand. If  slack  ex¬ 
ternal  flow  is  zero 
go  to  2.  Otherwise 
create  slack  arc  and 
store  data  in  the 
correct  position  in 
the  arc  list. 


READ , N, M:«0, SLACK: «N+1,XM:»0,  TOTAL :-0 


FOR  I:-l  TO  N 


V, 


B(I):-0 


NODE 


READ  I,BF,BS,CS 


I  -  0 


N| 


BF>0 


B(I):-BF 

XM:-XM+1,XK(XM):«I 

PARA(XM) :—  BF ,  TOTAL :  -TOTAL-BF 

L _ 

\ 

o 

• 

s 

BS 

>  o  yi 

J:-I,I:-SLACK, 

J: -5 LACK, LOWER :-0 

LOWER :-0, UPPER : -BS 

UPPER:— BS, 

COST :-CS, GAIN :-l 

COST :-CS, GAIN :-l 

(  ORIGSGKI.J, LOWER 

.UPPER,  COST,  GAIN)) 

- >NODE 


is  blank  go  Co  4. 
Else  store  Che  are 
in  arc  list. 

4. (EXT)  Place  each 
arc's  data  in  the 
correct  position 
on  the  arc  term¬ 
inal  list. 


ARC 


READ  1,J, LOWER, UPPER, COST, GAIN 

I  -  0  y 

5N: -SLACK 
I“N  :-0 

foRIGSGl (I, J, LOWER, UP PER, COST, GAIN 

- >ARC 

EXT 

LM:-M,M:-0 

/  / 

FOR  K:-l  TO  LM 

f  / 

u 

u 

u 

J:-T(K),M:-M+1 

Cterms(k.j)) 

RETURN 

Subroutine  READNDE(cont) 
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PURPOSE:  To  im¬ 
plement  Che  para- 
■eerie  node  flow 
augmentation  ap¬ 
proach  for  a  gener¬ 
alized  network . 

1. ( INITIAL)  Init¬ 
ialize  all  para¬ 
meters  to  zero. 

2.  (OPTIMAL)  If  the 
difference  between 
the  emulative  re¬ 
ference  nmber  and 
1  is  zero,  stop 
with  the  optimal 
solution.  Other¬ 
wise  go  to  3. 

3.  (ZERO)  Set  the 
y' s  to  zero.  Init¬ 
ialize  arc, node, 
and  cycle  counters 
to  zero. 

4.  (LEAVE)  For  each 
node  in  the  candi- 
node  list  deter¬ 
mine  if  its  para¬ 
metric  parameter 
is  less  than  zero. 
If  so  calculate  a 
reference  nmber. 
Store  the  smallest 
as  ME.  Set  IJ  to 
the  node  nmber 
and  TM  to  1  or  -1 
as  appropriate. 

Call  PATH  which 
will  create  LISA 
LISN,  and  RC  lists 
for  all  flow  aug¬ 
menting  trails. 

If  a  cycle(s)  has 
been  found  in  PATH 


ZERO 


Subroutine  PARANDE 


2 


co  update  Y's  for 
aodaa  on  che  ey¬ 
elet  a).  Call 
MFLOG  to  determine 
Cha  reference  nua- 
ber  and  che  are 
which  reaehea  cap¬ 
acity.  If  MF+ 
CUMRNO  >  1  reset 
MF  and  update  flows 
and  return.  If  not 
check  If  a  node 
will  be  infeasible 
at  next  iteration. 
If  so  update  flows 
and  return.  Other-* 
wise  update  flows. 

4 (ENTER)  Call  FIND- 
NDE  which  will  de¬ 
termine  the  enter¬ 
ing  arc.  If  no  arc 
can  be  found,  there 
is  no  feasible  sol¬ 
ution.  Return. 
Otherwise  go  to 
5. 

5. (Change)  Call 
PIV0T1G  to  change 
basis  and  up  date 
dual  variables. 


(  PATH(IJ.TM.XKE.XKK.PX.I1C.CK)T 


// 


_ 

Z 

s  0 

_ 

CYCLE G 

MF:*IE,ITER:-ITER+i 

C  MFLOG  ) 

MF  +  CUMRNO  >1.0 

MF  -  1 -CUMRNO 

P 

■ 

o 

/N 

(  FLOWG  ) 

(  FLOWG  ) 

MF  « 

» o  yk 

IIERD:- 

ITERD+1 

CUMNRO :■ 
CUMNRO+ 

MF 

C  FLOWG  ) 

RETURN 

RETURN 

- >ENTER 

ENTER 


CHANGE 
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